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Abstract 



Studying the acceleration and propagation mechanisms of Galactic cosmic rays can 
provide information regarding astrophysical sources, the properties of our Galaxy, 
and possible exotic sources such as dark matter. To understand cosmic ray accel- 
eration and propagation mechanisms, accurate measurements of different cosmic 
ray elements over a wide energy range are needed. The PAMELA experiment is a 
satellite-borne apparatus which allows different cosmic ray species to be identified 
over background. 

Measurements of the cosmic ray antiproton flux and the antiproton-to-proton 
flux ratio from 1.5 GeV to 180 GeV are presented in this thesis, employing the 
data collected between June 2006 and December 2008. Compared to previous ex- 
periments, PAMELA extends the energy range of antiproton measurements and 
provides significantly higher statistics. During about 800 days of data collection, 
PAMELA identified approximately 1300 antiprotons including 61 above 31.7 GeV. 
A dramatic improvement of statistics is evident since only 2 events above 30 GeV 
are reported by previous experiments. The derived antiproton flux and antiproton- 
to-proton flux ratio are consistent with previous measurements and generally con- 
sidered to be produced as secondary products when cosmic ray protons and helium 
nuclei interact with the interstellar medium. 

To constrain cosmic ray acceleration and propagation models, the antiproton 
data measured by PAMELA were further used together with the proton spectrum 
reported by PAMELA, as well as the B/C data provided by other experiments. Sta- 
tistical tools were interfaced with the cosmic ray propagation package GALPROP 
to perform the constraining analyses. 

Different diffusion models were studied. It was shown in this work that only 
current PAMELA data, i.e. the antiproton-to-proton ratio and the proton flux, arc 
not able to place strong constraints on propagation parameters. Diffusion models 
with a linear diffusion coefficient and modified diffusion models with a low energy 
dependence of the diffusion coefficient were studied in the % 2 study. Uncertainties 
on the parameters and the goodness of fit of each model were given. Some models 
are further studied using the Bayesian inference. Posterior means and errors of the 
parameters base on our prior knowledge on them were obtained in the Bayesian 
framework. This method also allowed us to understand the correlation between 
parameters and compare models. 
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Since the B/C ratio used in this analysis is from experiments other than PAMELA, 
future PAMELA secondary-to-primary ratios (B/C, 2 H/ 4 He and 3 He/ 4 He) can 
be used to avoid the data sets inconsistencies between different experiments and 
to minimize uncertainties on the solar modulation parameters. More robust and 
tighter constraints are expected. The statistical techniques have been demonstrated 
useful to constrain models and can be extended to other observations, e.g. elec- 
trons, positrons, gamma rays etc. Using these channels, exotic contributions from, 
for example, dark matter will be further investigated in future. 
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Introduction 



Outline of the thesis 

This thesis presents measurements of cosmic ray antiprotons performed with the 
PAMELA^ satellite experiment. Cosmic ray propagation models are studied by 
using the antiproton and proton data measured by PAMELA and measurements 
of the B/C ratio from other experiments. An overview of cosmic rays is given 
in chapter 1, including the acceleration and transport mechanisms of cosmic rays, 
detection techniques for cosmic rays, and knowledge we obtain from cosmic ray 
studies. Chapter 2 further details the possible processes during cosmic ray prop- 
agation in our Galaxy and summarizes the current status of previous studies on 
cosmic ray propagation. Chapter 3 describes the PAMELA experiment. The sci- 
entific objectives of the experiment are illustrated. The design and identification 
capabilities of all the sub-detectors are detailed. Chapter 4 identifies antiprotons 
from a large background of various cosmic ray species and reconstructs the antipro- 
ton flux and the antiproton-to-proton ratio by estimating the selection efficiencies 
as well as other correction factors. In chapter 5, by employing the antiproton and 
proton data from PAEMLA and the B/C ratio data from other experiments, the \ 2 
minimization method and the Bayesian inference are used to constrain cosmic ray 
propagation models. Finally, some discussion and outlook are given in chapter 6. 

The author's contribution 

My work on PAMELA started in September 2007 when I started my PhD position 
in the group of Particle and Astroparticle Physics at KTH. The first year as a 
PhD student was mainly focused on familiarization of the PAMELA experiment 
and the data analysis framework. A few months were spent on an analysis to 
study the separation capability of particles with equal charge but different mass in 
the calorimeter by using a p/3 (momentum-velocity) method based on the multiple 
scattering effect. The analysis was documented as a Collaboration note, but not 
described further in this thesis. 

a Payload for Antimatter Matter Exploration and Light-nuclei Astrophysics. 
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From the second year, I took part in the analysis of antiproton measurements. 
Building on the antiproton selection criteria developed by the PAMELA Collabora- 
tion, I started working on estimating the antiproton selection efficiencies and recon- 
structing antiproton flux. The efficiencies were derived by using different methods 
to fully understand the detector performance and possible systematic effects. This 
analysis also provided some fundamental information for the analysis of proton 
measurements, which is performed elsewhere in the PAMELA Collaboration. This 
work was presented in my licentiate thesis in April 2010, entitled "Measurements of 
cosmic ray antiprotons with PAMELA" . A part of this doctoral thesis concerning 
about the antiproton measurements, i.e. chapter 4, is selected from the licentiate 
thesis. 

After doing the data analysis on antiproton measurements, I focused on study- 
ing cosmic ray propagation models by using statistical methods. The source and 
propagation parameters charactering the injection primary cosmic ray spectrum 
and different propagation processes were constrained under the framework of dif- 
ferent propagation models. The GALPROP package which solves the transport 
equation numerically was used in my work to simulate cosmic ray propagation. 
I interfaced GALPROP with the statistical tools MINU1T and MULTINEST to 
perform a x 2 minimization analysis and a Bayesian analysis, respectively. The 
constraining capability of current antiproton and proton data from PAMELA data 
were demonstrated, as well as the upcoming PAMELA B/C data. Furthermore, in 
the Bayesian analysis I also produced the credible intervals on the parameters as 
well as on the predicted cosmic ray fluxes and flux ratios to understand the sta- 
tistical uncertainties on parameters and predicted fluxes as well as the correlation 
between parameters. 

My work has been presented at several PAMELA Collaboration meetings and 
international conferences and has been discussed in several publications. 

Publications 

• J. Wu, "Measurements of cosmic-ray antiprotons with PAMELA" , Royal In- 
stitute of Techology Licentiate Thesis (2010), ISBN: 978-91-7415-585-3. 

• J. Wu on behalf of the PAMELA collaboration, "Measurements of cosmic-ray 
antiprotons with PAMELA", Astrophys. Space Sci. Trans. 7 (2011) 225-228. 

• J. Wu et al., "Constraints on cosmic-ray propagation and acceleration models 
from recent data" , Proceedings of 32nd International Cosmic Ray Conference 
(2011). 

• O. Adriani et al., "PAMELA Results on the Cosmic-Ray Antiproton Flux 
from 60 MeV to 180 GcV in Kinetic Energy" , Physical Review Letters 105 
(2010) 121101. 
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• O. Adriani et al., "PAMELA Measurements of Cosmic-Ray Proton and He- 
lium Spectra", Science 332 (2011) 69-72. 

• O. Adriani et al., "Cosmic-Ray Electron Flux Measured by the PAMELA 
Experiment between 1 and 625 GeV", Physical Review Letters 106 (2011) 
201101. 

• O. Adriani et al., "The Discovery of Geomagnetically Trapped Cosmic-ray 
Antiprotons" , The Astrophysical Journal Letters, 737 (2011) L29. 

• O. Adriani et al., "Observations of the 2006 December 13 and 14 Solar Particle 
Events in the 80 McVn" 1 - 3 GcVn" 1 Range from Space with the PAMELA 
Detector", The Astrophysical Journal 742 (2011) 102. 

• O. Adriani et al., "A statistical procedure for the identification of positrons 
in the PAMELA experiment", Astroparticle Physics 34 (2010) 1-11. 

• O. Adriani et al., "Measurements of quasi-trapped electrons and positron 
fluxes with PAMELA", Journal of Geophysical Research 114 (2009) A12218. 

Presentations 

• 21st Nordic Conference in Particle Physics, Spatind, Norway. January 3-7, 

2010. Contributed talk, "The cosmic ray antiproton flux between 1 GeV and 
180 GeV measured by PAMELA". 

• A workshop on cosmic ray backgrounds for dark matter searches, Oscar Klein 
Center, Stockholm, Sweden. January 25-27, 2010. Contributed talk, "The 
cosmic-ray antiproton flux measured by PAMELA" . 

• 22nd European Cosmic Ray Symposium in Turku, Finland. August 3-6, 2010. 
Contributed talk, "Measurements of cosmic-ray antiprotons with PAMELA" . 

• 7th TeV Particle Astrophysics Conference, Stockholm, Sweden. August 1-5, 

2011. Contributed talk, "Constraints on cosmic-ray propagation and acceler- 
ation models from recent data" . 

• 32nd International Cosmic Ray Conference, Beijing, China. August 11-18, 
2011. Contributed talk, "Constraints on cosmic-ray propagation and acceler- 
ation models from recent data" . 
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Chapter 1 

Cosmic rays 



A general picture about the basics of cosmic rays is given in this chapter. In section 
|l.l| we review the characteristics of cosmic rays and discuss the questions concerning 
the hypotheses on cosmic ray production and propagation which remain unclear. 
Section |1.2| focuses on detection techniques for cosmic rays, from high altitude 
balloon-borne experiments or space missions to the ground based detection of ultra 
high energy particles. Finally, section [I~3| presents the use of cosmic rays as a tool to 
study aspects of astrophysics, dark matter, and the matter-antimatter asymmetry 
in the Universe. 



I. 1 Introduction to cosmic rays 

Cosmic rays are energetic charged particles from outer space that travel at nearly 
the speed of light and impinge on Earth from all directions. They are composed 
mainly of ionized nuclei, roughly 90% protons, 10% helium nuclei, and slightly 
under 1% heavier elements as well as electrons. In 1912 Victor Hess found that an 
electroscope discharged faster as he ascended in a balloon to altitudes up to 5 km [1] . 
He therefore concluded that cosmic rays arrived from outside our atmosphere and 
did not originate from decaying radioactive isotopes in the ground. Since their 
discovery, cosmic ray nuclei and electrons have been studied extensively [2] with a 
special emphasis on their main characteristics: the energy spectrum as well as their 
elementary composition and abundances. 

The overall energy spectrum of cosmic ray for energies > 10 10 eV, where solar 
effects are negligible, is well described by an inverse power law with evident features: 
the "knee" , at 4 x 10 15 eV, a not so evident second knee at ~ 10 17 eV, and a flatter 
supposedly extragalactic component at energies larger than about 5 x 10 18 eV (figure 

I I . 1 [ ) . Beyond this energy, where data become sparse, another steepening appears 
above 5 x 10 19 eV possibly due to the Greisen-Zatsepin-Kuzmin limit (GZK limit), 
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which is expected as cosmic ray protons with energies above 5 x 10 19 eV interact 
with cosmic microwave background photons to produce pions: 

a + P + 7T° 

p + JCMB^A -> n + n+ . 

The main features of the cosmic ray composition at low energies (< 10 14 eV) 
were known by 1950 and still remain unclear at higher energies. The relative 
abundances of cosmic rays are similar to the abundances of common elements in the 
Solar system, as shown in figure [TT2[ This consistency indicates that the composition 
of cosmic ray material injected into the interstellar medium (ISM) is very similar 
to that of the nebula that formed the Solar system. However, a striking difference 
can be seen between these two compositions. Chemical elements including Li, Be, 
B, F, Sc, Ti, V, Cr, Mn which are rare in the Solar system are many orders of 
magnitude more abundant in the cosmic rays. Since these elements are essentially 
absent as end products of stellar nucleosynthesis, they are generated as spallation 
products of abundant cosmic rays interacting with hydrogen or helium nuclei in the 
interstellar gas. For instance, Li, Be, B isotopes are mainly created from fragmented 
progenitors C, N and O nuclei. An example of the reaction contributing to boron 
production is 12 C + p — > 10 B + 3 He. 

Although a great amount of information on the composition and the energy spec- 
trum of cosmic rays on Earth has been gathered, fundamental questions concerning 
the origin of these particles, the mechanism through which they are accelerated to 
high energies, and the processes they undergoing before they arrive at the Earth 
remain unanswered. 



1.1.1 Cosmic ray sources and acceleration 

Since the 1960s, supernova remnants (SNRs) - the tattered, gaseous remains of 
supernovae - have been discussed as the breeding ground of Galactic cosmic rays 
for energies up to 10 15 eV. On average, about one supernova occurs in our Galaxy 
every 30 years, releasing 10 44 J in the form of kinetic energy in the ejecta. Therefore 
supernovae have enough power to energize the Galactic cosmic ray population at 
the observed level if there exists a mechanism for converting about 10% of the 
mechanical energy into relativistic particles. The knee in the cosmic ray energy 
spectrum presumably indicates a limit for the acceleration of cosmic ray protons 
by SNRs. It is argued that type II supernova surrounded by a dense stellar wind 
may be responsible for accelerating cosmic ray heavy nuclei up to energies about 
10 18 eV [5]. Above ~ 10 19 eV, the Galactic magnetic field would not be able to trap 
effectively even the heaviest elements of cosmic rays and an extra-galactic origin is 
required. Candidates such as external shocks in jets of active galactic nuclei and 
long gamma ray bursts have been proposed [5J [7J [5] as sources of ultra high energy 
cosmic rays. 

Once created in the sources, cosmic rays need to be accelerated and injected 
into the ISM. Diffusive shock acceleration (DSA) operating at expanding supernova 
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Figure 1.1. The energy spectra of cosmic rays (taken from [3])- Above 10 10 eV the 
spectrum shows a power-law behaviour. An obvious change in the slope is observed 
at the knee (4 X 10 15 eV) and at the ankle (5 X 10 18 eV). 
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Figure 1.2. The relative abundances of cosmic rays measured at Earth compared 
to the Solar system abundances (taken from [4], normalized to Si=100). 



shells is the most-favored mechanism for the production and acceleration of Galactic 
cosmic rays. The supernova remnants expand into the surrounding interstellar gas, 
compressing both the interstellar gas and magnetic field, producing a shock front. 
As the fast-moving charged particles move through the shocked gas, they diffuse 
by scattering on the contorted magnetic fields. Particles gain energy by bouncing 
between converging upstream and downstream regions around the shock front. This 
process naturally generates power-law energy spectra N (E) oc E~ a which is the 
striking characteristic of cosmic rays, with a = 2 for large sonic Mach numbers 
[HI EH EU HH IS]. In a more realistic picture, cosmic rays being accelerated can 
cause streaming instabilities and generate hydromagnetic waves which make the 
acceleration a non-linear process. A deviation from a = 2 can then occur. However, 
the injection spectrum remains poorly known since it does not only depend on the 
instantaneous spectrum of particles being accelerated at a shock, but also relates 
to how and when accelerated particles are released into the Galaxy, as well as the 
details of the interplay between accelerated particles, magnetic field amplification 
and shock dynamics. Calculations of non-linear DSA (NLDSA) models can either 
predict a hard injection spectrum with a spectral index less than ~ 2.1 — 2.15 
[T2J [13l [14] , or produce steeper spectrum up to a ~ 2.5 [15]. 
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Despite the appeal of the SNR conjecture, verification from observational evi- 
dence is needed. The problem is that cosmic rays are deflected and isotropized by 
the Galactic magnetic field and as a result the actual position of acceleration sites 
can not be extrapolated from their arrival direction. Thus some other tools are 
required to test the supernova paradigm. Significant progress has been achieved 
in recent years by keV X-ray and GeV to TeV gamma ray observations of young 
SNRs, providing very useful information about cosmic ray acceleration by super- 
nova shocks. Since the acceleration of cosmic rays in SNRs must be accompanied 
by copious gamma ray emission due to the decay of neutral pions produced in nu- 
clear collisions between relativistic nuclei and the background gas atoms, gamma 
ray detection is a good tracer of cosmic ray accelerators. Young SNRs which have 
strong shocks and can actively accelerate particles to the highest energies are usu- 
ally chosen to be targets to investigate acceleration processes. 

Many young SNRs exhibit shell-like morphologies at different wavelength bands. 



Examples of Tycho PU and RX J1713.7-3946 [171 ITS] are shown in figure[L3j While 
the non-thermal X-rays detected in the shells of SNRs are generated by electrons 
via synchrontron processes, the mechanism responsible for the gamma ray emis- 
sion is still under debate. Two scenarios have been proposed. Hadronic models 
connect gamma rays with neutral pion decay following proton-proton interactions 
while the leptonic models suggest gamma rays are generated through inverse Comp- 
ton scattering by the same populations of electrons interpreting the X-ray emission 
[19"1 [20l |2T] . Very recently, GeV gamma ray emission from RX J1713. 7-3946 was 
measured by Fermi-LAT [22] and disfavors the 7r°-decay mechanism. A very hard 
photon spectrum was observed which agrees well with the prediction of leptonic 
origin. However, the observed photon flux could still be reasonable in hadronic 
models considering low-density hot bubbles around the SNR shocks [33] or interac- 
tion between shocks and interstellar clouds [24l [25] . Another young SNR, Tycho, 
newly detected in GeV energies by Fermi-LAT [26] and in TeV energies by VERI- 
TAS [57], strongly supports the hadronic scenario, as shown in figure 1.4 However, 
both scenarios could be adapted to the experimental data under the assumption of 
a SNR environment with non- uniform magnetic fields |28j . Neutrino observations 
with km 3 -class detectors such as IceCube |3U] or KM3NeT [3T] may improve 
our confidence in the hadronic mechanism, since high energy neutrinos are created 
mainly in the decay of charged pion mesons produced in collisions of cosmic ray 
protons with nuclei in the ambient gas. 



1.1.2 The journey of cosmic rays from the source to the 
Earth 

No matter where and how the cosmic rays were produced and accelerated, they sub- 
sequently propagate through the Galaxy before reaching Earth. A common view 
derived from observations hints that cosmic rays travel in a confinement volume 
with an average residence time ~ 10 7 years. The amount of matter traversed by 
cosmic rays is estimated to be less than the density of the disk, indicating that 
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Figure 1.3. Left: Three-color composite image of Tycho's SNR observed by Chan- 
dra (taken from Hi : 0.95-1.26 keV emitted from Fe L-shell (red), 1.63-2.26 keV 
emitted from Si K-shell (green) and 4.1-6.1 keV continuum (blue). Right: RX 
J1713. 7-3946 as seen by HESS (colors) and by ASCA in the l-3keV energy band 
(contours) (taken from 1181 ). The image is smoothed with a Gaussian of 2 and the 
linear color scale is in units of excess counts per smoothing radius. 




10 8 10 9 10 10 10 11 10 12 10 

E [eV] 



Figure 1.4. The spectrum of gamma ray emission from Tycho's SNR measured by 
Fermi-LAT and VERITAS, compared with different theoretical contributions. Taken 
from [32]. 
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cosmic rays are trapped mostly in low-density regions. Several complex phenom- 
ena might occur during propagation. It is believed that during their residence 
time in the Galaxy, cosmic rays diffuse randomly due to the irregularities in the 
Galactic magnetic field. The diffusion process was proposed to explain cosmic ray 
confinement in the Galaxy and the observed isotropy. The details of the Galactic 
magnetic field structure is not well understood. Assuming the average magnetic 
field strength is B, a particle scattering on magnetic field irregularities with weak 
random fluctuations SB « B can be treated in the quasi-linear theory of plasma 
turbulence. Theoretically, different types of spectral energy density of interstellar 
turbulence have been proposed. The favored ones are Kolmogorov-type |33j and 
Kraichnan-type [34j spectra. However, current data do not allow us to distinguish 
between these different turbulence types. 

In addition to diffusion, other processes could also play a role in cosmic ray 
transport. Cosmic rays could possibly be convected if the medium responsible for 
diffusion is moving away from the disc, i.e. a galactic wind is present. The scat- 
tering of cosmic ray particles on magnetized plasmas in the ISM causing stochastic 
acceleration could also happen, but cannot serve as the main mechanism of cosmic 
ray acceleration. Moreover, when charged cosmic ray nuclei travel in the ISM, they 
undergo nuclear destruction due to fragmentation and unstable nuclei decay to sta- 
ble nuclei. Through these processes, secondary cosmic rays are created as spallation 
products of primary progenitors. Additionally, energy losses arise from interactions 
such as ionization and Compton scattering, which dominate for cosmic ray nuclei. 
Cosmic ray electrons, however, do not only lose energy by virtue of interactions 
with the ISM but also with the Galactic magnetic field or the interstellar radiation 
field. A more detailed description of the propagation processes in the Galaxy will 
be a focus of Chapter [2] 

Finally, before arriving at Earth, cosmic rays are affected by the outstreaming 
particles ejected from Sun and the geomagnetic field. The Sun emits low energy 
particles in the form of a fully ionized plasma called the solar wind, dominating in a 
cavity known as heliosphere as shown in figure [T75| The solar wind has a supersonic 
speed of about 400-800 km/s, flows outward and decreases to subsonic flow at the 
termination shock. Beyond this, the solar wind which carries the spiraling inter- 
planetary magnetic field is turned toward the heliotail. At larger radial distances, a 
surface called the heliopause is reached, separating the solar material and the solar 
magnetic fields from the interstellar material and the interstellar magnetic fields. 
Interstellar ions are diverted around the heliosphere. An outward pointing bow 
shock may also be formed beyond the heliosphere. 

The solar wind prevents low energy cosmic rays from penetrating the heliosphere 
and modulates the cosmic ray energy spectra. This phenomenon is called solar 
modulation and was developed originally by Parker |36| . varying according to the 
11 year solar cycle. The greater the solar activity, fewer cosmic ray particles can 
get into the heliosphere, as shown in figure [T~6} The solar modulation is determined 
by four mechanisms, including convection by the outward solar wind flow, diffusion 
in a turbulent heliospheric magnetic field (HMF) carried by the wind, drift due 
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Figure 1.5. A schematic diagram of the heliosphere. Taken from |35| . 



to the gradients, curvature and current sheet of the HMF and adiabatic energy 
changes. A simple but frequently used model is the force-field approximation |37j 
which depends on a single parameter, the modulation potential For a nucleus 
with charge Z, mass m and atomic number A, its interstellar flux Jjg is modulated 
to the top-of-atmosphere flux Jtoa by the relation 

Jtoa(e) = {E+ %^y_ m2 J IS (E + W*), (i.i) 

where E is the total energy of the nucleus. The modulation potential $ is deter- 
mined by fitting the observed spectrum above the atmosphere with the assumed 
interstellar spectrum. Although the force-field approximation is useful in most 
cases, it is worth to note that this approximation cannot consider any charge-sign 
dependence of solar modulation indicated in experimental data |38l 1183) . Drift 
models are suggested in literature [301 S3 HH S3] which can produce a clear charge- 
sign-dependent modulation. For example, during the A < polarity cycles, i.e. 
when the HMF is directed toward the Sun in the northern hemisphere, the neg- 
atively charged particles will drift inward primarily through the polar regions of 
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the heliosphere and the positively charged particles will drift primarily through the 
equatorial regions of the heliosphere. Nevertheless, the realistic time-dependent 
modulation could be very complex and needs further investigation [43] . Cosmic 
ray nuclei with energies larger than about 10 GeV/n are not sensitive to the solar 
activity [33] . 
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Figure 1.6. Variation of cosmic ray neutron intensity and the solar activity rep- 
resented by the sunspot numbers. High cosmic ray intensity corresponds to low 
sunspot activity, and vice versa. Taken from |45j . 



After penetrating the heliosphere, low energy charged particles are deflected by 
the geomagnetic field, which is the last obstacle for cosmic rays on their way to 
the top of Earth's atmosphere. A charged particle traverses this magnetic field in 
a curved path and a minimum rigidity (momentum per unit charge) referred to 
as the cutoff is required to penetrate the geomagnetic field. The cutoff rigidity, 
varying with the geomagnetic position and the approaching direction of cosmic ray 
particles, was first treated by Stoermer who approximated the Earth's magnetic 
field as a dipolar field [46]. For particles incident vertically towards the center of 
the magnetic dipole, the Stoermer vertical cutoff (SVC) can be written as [47] 



p> 14.9Zcos 4 A GeV/c, 



(1.2) 



where A is the geomagnetic latitude. Consequently, the detected cosmic ray inten- 
sity will be lower at the magnetic equator and higher at the magnetic pole as the 
geomagnetic cutoff value is largest at the equator and diminishes closer to the poles. 
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The SVC is often the reference quantity calculated and is used as an effective aver- 
age over all arrival directions. However, the SVC has limited accuracy because the 
realistic geomagnetic field does not obey an ideal dipole geometry but is offset by 
some 400 km from Earth's center and has higher order components. Furthermore, 
Stoermer's theory allows the penetration of charged particles with trajectories that 
would go through Earth and generally underestimates the cutoffs. This problem is 
called the Earth's shadow and was first addressed by Vallarta [48], who showed that 
a range of magnetic rigidities exist above the Stoermer cutoff where the penumbral 
shadow of Earth casts a broken pattern of allowed or forbidden bands of magnetic 
rigidity. More reliable and precise determination of the geomagnetic field can be 
done by tracing trajectories of cosmic rays in higher order geomagnetic field models 
091 [50]. 



1.2 Detection techniques 

In order to understand the nature of cosmic rays, many experiments have been 
performed during the last century, producing a large amount of observational data. 
Different kinds of detectors arc used to detect cosmic rays depending on the energy 
of interest. Direct detection experiments record cosmic rays directly, while indi- 
rect ones measure the secondary showers initiated from the incident cosmic rays 
interacting with atmosphere. 

Direct detection is used to study particles below 10 15 eV for which the flux of 
particles is sufficiently large that individual primary nuclei can be studied by instru- 
ments carried in high-altitude balloons or in space. The purpose of direct detection 
is to discriminate the incoming cosmic ray particles and to measure their abun- 
dances and energies. Various types of detector are utilized, such as magnetic spec- 
trometers, calorimeters, transition radiation detectors, scintillators or solid state 
detectors, Cherenkov counters and time-of-flight systems. A number of these de- 
tectors are appropriately assembled as a package either in high-attitude balloon 
experiments such as MASS91 [51], CAPRICE E21[53], TRACER [54], ATIC [55], 
BESS [H] and CREAM [57 , or in space-based experiments, for example Spacelab 2 
[SS], HEA03-C2 |52|, ACE-CRIS |gg, AMS [SUES] and PAMELA [53]. In general, 
balloon-borne experiments allow multiple flights with a moderate budget and can 
provide a prototype test which can be further employed in space. However, the 
exposure time they can provide is up to 42 days (CREAM-1 performed in 2004 
|57|). which is restricted mostly by the wind and limited resources on-board. Space 
missions are more expensive and risky, but highly increase the statistics benefiting 
from much longer exposure time and reduce the systematic uncertainties caused 
by the interference of cosmic rays with the residual atmosphere above balloons. A 
part of this thesis focuses on the data analysis of the satellite-borne experiment 
PAMELA, which will be described in more detail in chapter [3j 
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Very high energy (above ~10 14 eV) cosmic rays are extremely rare, for example 
the flux is only 1 km _2 sr~ 1 year _1 above 10 19 eV (see figure 1.1), and only ground- 
based experiments with huge effective areas and long exposure times can hope 
to acquire a significant statistical sample. The ground-based experiments exploit 
the atmosphere as a giant calorimeter. An incident cosmic ray particle interacts 
with air molecules, mainly oxygen and nitrogen, and produces a cascade of lighter 
particles, spreading out over large areas, called an extensive air shower. Rather 
than detecting the primary cosmic rays directly, ground-based detectors detect 
the remnants of the atmospheric cascades of particles initiated by the primary 
particle. Composition and energy information of incident particle species can be 
derived from the EAS properties based on hadronic models. Several techniques are 
used in current instruments, ranging from direct sampling of secondary particles 
in the shower to measurements of fluorescence from atmospheric nitrogen excited 
by the charged particles, and radio emission emanating from the air shower. Some 
experiments employing one or more of these techniques, are AGASA [64j . HiRcs 
E5|, Auger [S3, KASCADE jUj and TA JSSJ. 



1.3 Cosmic rays as observational tools 

From the various elemental cosmic ray data, measured by different experiments, 
the principle astrophysical issues concerning cosmic ray acceleration and propaga- 
tion mechanisms can be investigated. Moreover, cosmic ray observations provide 
potential to help us understand topics such as the nature of dark matter and the 
apparent matter-antimatter asymmetry in the Universe. 



1.3.1 Astrophysics 

The energy spectra of cosmic rays, extending over a wide energy range from tens 
of MeV/n to EeV/n, provide a useful means to probe the properties of cosmic 
ray sources, acceleration mechanisms, propagation processes in the Galactic halo 
and the interstellar environment itself. At energies higher than tens of GeV, the 
observed abundances are affected by the injection spectrum from the sources, the 
diffusion in the Galactic magnetic field and the nuclear interactions in the Galaxy. 
The low energy tail, however, also has a contribution from other phenomena, such 
as convection, reacceleration and heliospheric physics. Therefore, tracing back from 
cosmic rays observed at Earth, we can effectively investigate the processes happen- 
ing before cosmic rays reach Earth. In addition, the properties of the ISM and the 
structure of the Galactic magnetic field can then be better understood. 

As indicated in figure [OJ elements like Li, Be, B are secondary nuclei produced 
by primary cosmic rays interacting with the interstellar gas. Therefore, the rela- 
tive abundances of secondary nuclei shed light on the properties of matter in the 
Galaxy. Measurements of secondary-to-primary ratios are useful probes of cosmic 
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ray transport since they mainly depend on the mean amount of interstellar mat- 
ter that primaries have encountered before reaching Earth rather than the source 
spectrum of the progenitors. The B/C ratio has been considered as one of the most 
important quantities for decades, as B is entirely secondary and its main progen- 
itors C and O are primaries directly produced in the SNR. The B/C ratio is also 
the best measured secondary-to-primary ratio since it depends on the elemental 
separation capability of the detector but not on the isotopic separation capability 
which are important for the ratios like 2 H/ 4 He. Figure [L7] shows the measured 
B/C ratio compared to some models, which cannot be distinguished using the data 
listed here alone. Apart from B/C data, other quantities such as 2 H/ 4 He, 3 He/ 4 He 
and p/p are also secondary-to-primary ratios which are useful to probe cosmic ray 
transport processes [175] . 
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Figure 1.7. Measured boron-to-carbon ratio (B/C) compared with four models, 
which are modulated with the solar modulation potential <1> = 500 MV. Taken from 

The other major constraint on propagation models comes from radioactive 
species, which are unstable nuclei undergoing radioactive decays such as /3 decay 
and electron capture. Especially useful species are secondary radioactive isotopes 
for which no extra contribution from sources need to be accounted for. The ratios 
of unstable to stable isotopes of secondary nuclei tell us the global properties of the 
Galaxy through the surviving fraction of unstable isotopes in the Galaxy. A com- 
bination of secondary-to-primary ratios and radioactive isotope ratios allows one to 
derive the size of Galaxy halo. The most notable unstable nucleus is 10 Be, which 
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is best measured and has a lifetime of ~ 3.9 x 10 6 years for /? decay, comparable 
with the escape time of cosmic rays in the Galaxy. Other long-lived radioactive 
nuclei such as 14 C, 26 Al, 36 C1 and 54 Mn also provide constraints on cosmic ray 
propagation [70II7T] . 

The spectra observed on Earth are affected by a combination of acceleration 
and propagation. Generally, while the propagation processes can be understood 
using secondary cosmic rays, the information on the acceleration can be derived 
from the primary cosmic ray spectra |72j . Propagation and source parameters are 
degenerate. Simultaneously fitting secondary-to-primary ratios as well as primary 
fluxes allows us to explore both source and propagation mechanisms |73j . 

1.3.2 Dark Matter 

The existence of dark matter (DM) is motivated by a wealth of observational evi- 
dence, including galactic rotation curves, gravitational lensing, the anisotropies of 
the cosmic microwave background, and primordial light element abundances. The 
approximate distribution of DM, which constitutes about a quarter of the mass of 
the Universe, can be deduced from its gravitational effects, but its nature and mi- 
crophysical properties remain one of the great unsolved problems of physics |74U75) . 
The lack of observation of DM particles indicate that DM particles are primarily 
non-baryonic which only interact through the weak force and gravity. In addition, 
cold DM is necessary to explain structure formation in Universe, since relativistic 
DM moves too quickly to clump together on small scale of galaxies. One of the 
most common proposed candidates is referred to as a Weakly Interacting Massive 
Particle (WIMP) [73 [77]. DM signals are possible to be detected directly by cre- 
ating DM particles in accelerators or searching for the scattering of DM particles 
off atomic nuclei within a detector, and indirectly by gathering information from 
WIMP annihilation products. 

Indirect detection of DM is based on the search for anomalous features in energy 
spectrum in cosmic rays due to WIMP annihilation in the Galactic halo, on the top 
of the expectation from standard astrophysics background. Since matter dominates 
cosmic rays, antimatter (p, D, e + ), gamma ray and neutrino channels have better 
potential to probe dark matter. The spectral distortion of these components over 
the astrophysical background may give evidence for dark matter. 

The DM signal prediction depends on the models in which the properties of 
DM particles and their interaction strength with Standard Model (SM) states are 
assumed. There are numerous models with a WIMP DM candidate. One popu- 
lar framework is the Minimal Supersymmetric extension to the Standard Model 
(MSSM), in which the lightest supersymmetry particle, known as "neutralino" , is 
stable and therefore provides a good candidate. Another widely discussed scenario 
is universal extra dimensions introducing a tower of Kaluza-Klein partners for ev- 
ery SM particle in which the lightest Kaluza-Klein particle is stable and considered 
as a good DM candidate. For any particular model, the annihilation products of 
DM candidates can be predicted. A robust and accurate estimation of the cosmic 
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ray background contribution as well as the propagation of DM annihilation prod- 
ucts is required to clarify DM models and consequently help us to understand the 
properties of DM particles. 

Recently, the positron fraction first reported by PAMELA [78] and then con- 
firmed by Fermi-LAT [79] show an unexpected excess above 10 GeV over the predic- 
tion of propagation models, triggering many theoretical interpretations including 
a contribution from DM (see, e.g. [SOI IHH H21 IB3 IHH IBS]). The antiproton flux 
is also an interesting means for indirect DM detection as they are inevitably pro- 
duced whenever it is kinematically possible and the final states of DM annihilation 
contain quarks or gauge bosons. One example is shown in figure |1.8[ which shows 
that the antiproton flux at energies of tens of GeV resulting from annihilation of 
high-mass neutralinos could be more than an order of magnitude above the flux of 
secondary antiprotons and thus could be observed. The PAMELA measurements 
of antiproton flux presented in this thesis, which extends to an energy of about 
180 GeV provides an important test for these models. 

1.3.3 Matter-antimatter asymmetry 

The matter-antimatter asymmetry in the Universe is one of the most important 
and puzzling questions indicated from cosmic rays observation. So far, no sizable 
amounts of antimatter has been observed. The asymmetry is also inferred from the 
baryon-to-photon ratio in the cosmic microwave background radiation, i.e. t\b ~ 
10~ 9 [87 . While fundamental theories for elementary particles predict the same 
laws for matter and antimatter, CP-violation and baryon number non-conservation 
have been proposed to explain the depletion of antimatter |88j . However, particle 
experiments do not support large levels of violations |89j . Domains of antimatter 
in our Universe are suggested in (90j [91] [92] ■ 

The detection of antinuclei with charge \Z\ >= 2 would constitute a smoking 
gun if they can be found in future experiments since the secondary production of 
antinuclei is negligible due to the extremely small production probability in the ISM 
through spallation (e.g. 10 -13 for 3 He). Inferred from cosmic ray matter compo- 
sition, cosmic ray antihclium are the most possible detectable antinuclei compared 
to other species. So far, no He has been detected and only an upper limit on He/He 
has been reported by experiments. Cosmic ray antiprotons and positrons, which 
are measured with much higher statistics than antihelium nuclei, could also provide 
signals on primordial antimatter sources. However, studies could become compli- 
cated since contributions from such as non-standard astrophysical sources and dark 
matter may also give an excess on the cosmic ray antiproton spectrum or positron 
spectrum. 



In summary, this chapter reviews the fundamental issues of cosmic rays. A lot 
of questions about their nature, origin and propagation are still unanswered. Ac- 
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Figure 1.8. A primary antiproton flux assuming the annihilation from neutralino 
from MSSM with a mass of 964 GeV (dotted line) compared with experiment results. 
The solid lines show the upper and lower limit of calculated flux of interstellar sec- 
ondary antiprotons by Simon el al. The dashed line shows the theoretical calculation 
of interstellar secondary antiprotons by Bergstrom & Ullio. All the references for 
the experiment results and theoretical models can be found in |86| . 

curate measurements of individual cosmic ray elements are necessary to study the 
cosmic ray acceleration and propagation phenomena. Additionally, dark matter 
and matter-antimatter asymmetry can be inferred from cosmic antimatter mea- 
surements. As outlined in this chapter, the study of cosmic ray propagation plays a 
key role in understanding the processes occurring in our Galaxy. The propagation 
processes briefly described here will be elaborated in chapter [2j 
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Chapter 2 

Cosmic ray propagation 



This chapter discusses general questions related to the propagation of cosmic rays 
with energies up to 10 15 eV. It starts with an overview of the basics of cosmic ray 
propagation, such as energy losses/gains and nuclear interactions of cosmic rays 
with the interstellar medium. In section 2.2 the transport equation is constructed 
in the form of the continuity equation and Fick's law. All relevant parameters used 
to characterize the cosmic ray propagation processes are described in section |2.3| 
Section [2~4] discusses different approaches to solve the transport equation. The final 
section summarizes the current status of studies on cosmic ray propagation. 



2.1 Basics of cosmic ray propagation 

As mentioned in chapter [I] due to insufficient information on the properties of the 
ISM and on the structure of the Galactic magnetic field, the specific mechanisms 
of cosmic ray propagation are not known yet. All our knowledge is developed and 
constructed semi-empirically based on cosmic ray observational results. So far it 
is generally supposed that the diffusion process with possible reacceleration and 
convection can be crucial during cosmic ray propagation. 

2.1.1 Diffusion 

From cosmic ray observations, especially secondary-to-primary ratios (e.g., B/C, 
p/p) and the unstable-to-stable isotope ratios of secondary nuclei (e.g., 10 Be/ 9 Be, 
26 A1/ 27 A1), the mean amount of matter (referred to as grammage X) traversed by 
cosmic rays and their escape time T esc from our Galaxy can be established. The 
grammage X is found to be about 5 g/cm 2 and the escape time r esc is estimated 
to be tens of million years. This suggests that cosmic rays travel in a confinement 
volume with average gas density p gas about 0.3 protons/cm 3 , deduced from the 
relation X — J vp gas T esc where v is the particle velocity. Since the Galactic plane 
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has an average gas density about 1 proton/cm 3 , it appears that cosmic rays must 
spend most of their time in low-density regions of the ISM, which could be either 
a hot coronal phase of the ISM and/or refer to a Galactic halo surrounding the 
disk with low gas density. Radio measurements support the halo hypothesis since 
significant amounts of synchrotron radiation are detected far away from the galactic 
disk and may be emitted from cosmic ray electrons (see [1] and references therein) . 

A process which can confine cosmic rays inside the Galaxy and can send them 
back to the disk from the halo is needed. A natural hypothesis is that cosmic 
rays scatter on turbulence in the magnetic field. Cosmic rays form a plasma of 
ionized particles. On the microscopic level, the cosmic ray particles interact with 
the magnetohydrodynamic (MHD) waves arising in magnetized plasmas. If the 
interaction is resonant, particles are scattered by the waves leading to diffusion. 
The diffusion process can explain not only the cosmic rays' long travel time but 
also their highly isotropic distribution in the Galaxy. If there was no scattering, 
due to the particular location of the Solar system there should be more cosmic rays 
from sources towards the Galactic center and we would expect a strong anisotropy 
towards this direction. Such an anisotropy is not seen for cosmic rays with energies 
less than 10 15 eV which apparently is destroyed by multiple scattering of the cosmic 
rays on their path from the sources to us. 

Locally, cosmic ray diffusion occurs along the magnetic field lines and thus can 
be quite anisotropic. However, on scales larger than 100 pc, the wandering of cosmic 
rays on irregular magnetic field can make the diffusion isotropic and randomize the 
trajectories of particles. 

2.1.2 Energy losses and gamma ray production 

During propagation, cosmic ray nuclei and electrons interact with other constituents 
of the Galaxy and continuously lose energy. Meanwhile, electromagnetic radiation 
is produced in some interactions. By observing the emissions from radio to gamma 
ray frequencies, it is possible to probe cosmic ray propagation. 

For relativistic cosmic ray electrons, considering the content and structure of 
the Galaxy, the following interactions may occur: 

• ionization of neutral interstellar matter. 

• Coulomb scattering of individual plasma electrons in the fully ionized plasma. 

• Bremsstrahlung in the neutral and ionized medium. An electron is deflected 
in the electrostatic potential of an atom, ion or molecule, losing energy by 
emitting a 7-ray photon. 

• Inverse Compton scattering of the interstellar radiation field. In this interac- 
tion the target photon is scattered to higher frequencies by receiving part of 
the kinetic energy transferred from the relativistic cosmic ray electron. 
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• Synchrotron radiation in magnetic fields. Electromagnetic radiation is emit- 
ted when a charged relativistic particle travels in a magnetic held that is 
uniform on scales much larger than the gyroradius of the particle. 

For cosmic ray nucleons, the cross sections for electromagnetic interactions of 
cosmic ray nuclei are much smaller than those of electrons. Therefore all the elec- 
tromagnetic processes responsible for electron energy losses can be neglected. In- 
teractions of cosmic ray nuclei with cosmic photons only becomes important at 
energies > 10 17 eV and will not be considered here. The remaining contributions 
are interactions with ISM: 

• Ionization of atoms and molecules in the ISM. 

• Coulomb interactions with the ionized plasma. 



When these interactions occur, the original cosmic ray energy spectrum and 
propagation processes are affected. The contributions of different interactions to 
the energy loss are energy dependent. The energy loss timescales, which are as- 
sociated to the inverse of energy loss rates, are shown in figure |2.1| For nucleons 
and low energy electrons, the most important processes responsible for the energy 
loss are Coulomb scattering and ionization. But for electrons with energies higher 
than around 1 GeV, synchrotron losses become dominant. A complete summary of 
energy losses can be found in [53"] . 
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Figure 2.1. Energy-loss timescales of nucleons (left) and electrons (right) in neutral 
and ionized hydrogen. In the left figure, solid lines show ionization losses and dashed 
lines show Coulomb losses. In the right figure, BO (BI) means the Bremsstrahlung 
losses in the neutral gas (ionized gas). The curves are calculated based on the 
assumption of equal neutral and ionized gas number densities (0.01 cm -3 ), and 
equal energy densities of photons and magnetic field (1 eV cm -3 ). Taken from |93| . 
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2.1.3 Nuclear interactions 

When cosmic rays traverse the ISM, they may interact with an interstellar hydrogen 
or helium nucleus and initiate nuclear reactions. For a specific type of cosmic ray 
nucleus, several kinds of nuclear interactions can be discussed: 

• Inelastic scattering of cosmic ray nuclei with ISM atoms and molecules which 
results in destruction of the given species governed by the total reaction cross 
section of that species. 

• Spallation determined by the formation rate from each parent element. 

• Radioactive decay of unstable cosmic ray nucleons. 

• Radioactive spallation generated via the decay of parent unstable isotopes. 



2.2 The transport equation 

Taking into account the decisive role played by diffusion as well as other possible 
interactions, the cosmic ray transport equation can be built by incorporating the 
continuity equation with Fick's law |94j . 

The fundamental continuity equation can be written as: 

ON 

^— V-J + ft (2.1) 

where N is the number density, J is its current generated due to a spatial gradient in 
the density N and q is the source term. Assuming that the diffusing particles obey 
Fick's law, J = —D\7N, where D is the diffusion tensor, the continuity equation 
leads to the diffusion approximation: 

°-g = V • (DVn) + q. (2.2) 

A rather general equation for cosmic ray species i can be constructed by taking 
into account all the relevant processes in addition to diffusion: 

• Continuous energy losses 

f = ->".), M 

where bi = dE/dt is the first order of energy loss; 

• Nuclear destruction 

= -nvo-iNi, (2.4) 

where n is the density of the interstellar gas, v is the particle velocity and cr^ 
is the inelastic scattering cross section of a nucleus of type i with nuclei of 
the interstellar gas; 
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Spallation from heavier nuclei 

m 

dt 



J2 nvo-ijNj, (2.5) 



where 5ij is the production cross section of nuclei of type i from heavier nuclei 
of type j; 

• Radioactive decay 

where Ti is the lifetime of a nucleus of type i; 

• Radioactive spallation 

dt 

rrij >rrii J 

where Ty is the lifetime of a nucleus of type j decaying radioactively to a 
nucleus of type i; 

Adding up all these terms, the diffusion equation can be written as: 
dNi 



ON, 1 , s 



77i ^ >m^ 



(2.8) 



2.2.1 Transport of cosmic rays by convection 



While diffusion is necessary to explain the high degree of isotropy and confinement 
in the Galaxy, other processes may also be of importance. Particularly, it is very 
likely that in our Galaxy there is large-scale motion of the interstellar gas with a 
"frozen" magnetic field, caused by the stellar activity and the energetic phenomena 
associated with the late stage of stellar evolution. This is referred to as convection 
or galactic wind. Cosmic rays are carried by the wind "as a whole" with some 
velocity V c outwards from the Galactic plane. Galactic winds are found in many 
galaxies |95j . It is natural to propose that supernovae also power a similar wind in 
the Milky Way [Ml M M M UM ■ 

Observational support comes from Galactic diffuse soft X-ray emission mea- 
sured by ROSAT, which can be interpreted by assuming the presence of a strong 
galactic wind in our Galaxy [101) . Convection adds a term - ( V • V c ) Ni to the 



diffusion equation (equation 2.8) and causes adiabatic energy losses, of the form 
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2.2.2 Transport of cosmic rays by reacceleration 

Though cosmic rays lose energy during propagation in the Galactic environment, 
they may also gain energy via stochastic acceleration. Above a few GeV/n, the frac- 
tion of secondary nuclei decrease as energy increases, which indicates that higher 
energy cosmic rays traverse less amount of matter (i.e. spend shorter time) in the 
Galaxy than lower energy ones. If acceleration only occurred together with frag- 
mentation, the expected time spent to accelerate cosmic rays to higher energies 
would be longer. Hence, cosmic rays are mainly accelerated before their propaga- 
tion, as discussed in chapter [JJ where diffusive shock acceleration at the shock- wave 
fronts in SNRs are proposed as the main mechanism of cosmic ray acceleration in 
the Galaxy. However, this does not exclude the possibility that cosmic rays might 
experience some additional acceleration after being injected from sources. Due to 
relativistic cosmic rays scattering on magnetic turbulence in the interstellar hydro- 
dynamical plasma, some weak stochastic acceleration is almost unavoidable. This 
can be called reacceleration. 

Reacceleration may be significant at low energies to explain the peaks of B/C 
ratio around 1 GeV/n as shown in figure [L~7{ but should only slightly distort the 
ratios above few GeV. Reacceleration leads to a second order energy gain, which 
adds a term ^fi 2 D pp ^^r to the transport equation, where D pp is the diffusion 
coefficient in momentum space. 

2.2.3 A full transport equation 

A schematic view of cosmic ray transport including all the most important propa- 
gation steps is illustrated in figure [272] By adding the convection and reacceleration 
terms, the full transport equation can be written as: 



m 

dt 



£ 

rrij >irii 

l 

<7j H 

Ti 



N j + V 



DVN % - V a Ni 



N l + ^D pp ^-^ bl N,-^ P - Nl 
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(2.9) 
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Figure 2.2. Schematic view of the propagation of cosmic rays in our Galaxy. After 
accelerated by SNR shock waves, cosmic rays suffer a combination of propagation 
processes in our Galaxy, including diffusion, reacceleration, convection, nuclear in- 
teraction, radioactive decay and energy losses. Taken from |116| . 

2.3 Parameter description 

In the framework of the diffusion equation, it is possible to use observational data 
to study the properties of cosmic ray composition, abundances, anisotropy and to 
determine the composition at the sources. By combining numerous experimental 
facts within certain models, the propagation parameters which will be further de- 
scribed in this section can be investigated. This allow us to better understand the 
related transport processes. 

Spatial diffusion 

Diffusion is a result of cosmic ray particles interacting with MHD waves. While our 
knowledge on the structure of the Galactic magnetic field is limited, we generally 
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assume that diffusion mainly takes place along the mean magnetic field direction. 
Charged cosmic rays scatter mainly on resonant magnetic field fluctuations and 
the diffusion coefficient is estimated to follow a rigidity power law according to the 
quasi-linear theory |103l 169] . The diffusion coefficient is usually represented by D xx 
and assumed to have the form: 

D xx = D p(-£\ , (2.10) 

where Dq is the normalization at reference rigidity po, linked to the fluctuation 
level of the hydromagnetic turbulence; the factor ft = v/c is the particle velocity 
and 5 the spectral index of diffusion coefficient related to the spectral index of 
turbulence spectrum. The rigidity p is usually used as the kinematic variable instead 
of momentum p. The free parameters concerning diffusion are D and 5. 



Reacceler at ion 

The energy gain through reacceleration is a result of diffusion in momentum space. 
The associated diffusion coefficient in momentum space D pp is taken from the model 
of minimal reacceleration by interstellar turbulence and is correlated to the velocity 
of disturbances in the hydrodynamical plasma, called the Alfven velocity. D pp is 
related to the spatial diffusion coefficient D xx with the expression [104] : 



Av 2 A p 2 
3(5 (4 - S 2 ) (4 - 6) D x 



where va is the Alfven velocity - the main free parameter related to reacceleration. 



Convection 

Considering the convection mechanism, through which cosmic rays can be trans- 
ported in bulk away from the Galactic plane, the convection velocity V c (z) is the 
quantity used to describe the convective wind. It is usually assumed that the ve- 
locity varies linearly with the distance from the Galactic plane z = as 

Vc (z)=V(0) + ^z, (2.12) 

in which 1^(0) is usually taken to be zero for a simplicity and dV/dz is the main 
free parameter. Besides, some studies |105|. 1116] assume a constant velocity in 
order to make the transport equation analytically solvable. Nevertheless, detailed 
information on the convection velocity is still unknown. 



Source term 

As well as transport processes, the source term is indispensable in order to describe 
cosmic ray data. In section [Li] it was stated that SNRs arc believed to be the main 
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sources of primary nuclei. For a cosmic ray species the injected density is assumed 
to be a power law in momentum p (or rigidity p) as expected from diffusive shock 
acceleration theory: 

q l {p) oc p~ v oc p~ v . (2-13) 

The general form of the average source term depends not only on the point- 
source injection spectrum but also on the spatial distribution of sources f{R, z) in 
the Galaxy: 

q i (p,f) = N i f(R,z)p- v , (2.14) 

where 2Vj is the normalisation abundance for the cosmic ray species i. The free 
parameters related to the source terms are the normalisation abundance iVj and 
the injection index v. 

The SNR distribution in the Galaxy is very poorly determined by radio surveys 
due to the small sample available and selection effects |106L 1107] . Pulsars could 
be a useful tracer of the SNR distribution since they are born in the core collapse 
of supernovae. Large samples of pulsars can be obtained, but could be biased by 
distance and interstellar dispersion uncertainties |108j . Moreover, the distributions 
of SNRs and pulsars as a function of galactocentric radius are both steeper than 
the distribution of cosmic ray sources chosen to reproduce the EGRET 7-ray data 
[55] , An enhancement of molecular gas in the outer Galaxy was proposed in |109| 
to moderate this 7-ray gradient problem but is disfavored by Fermi-LAT data [110) . 

As suggested in |lllj , the radial dependence of the SNR distribution can have 
the form: 

where Rq is the radial distance of the Sun from the Galactic center. Different values 
of the parameters (a, fj) have been adopted in the literature. For instance, (1.69, 
3.33) was found in [TUS] and (2.00, 3.53) in [W] to model the SNR distribution. The 
combination (2.35, 0.654) was obtained in |108) to fit the pulsar distribution whereas 
(0.5, 1.0) was determined to reproduce observed EGRET 7-ray based gradient [112) . 
These distributions are shown in figure |2~3| More sets of favored values can be found 
in [113] and references therein. 
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R (kpc) 

Figure 2.3. Cosmic ray source density as function of galactoccntric radius R, 
normalized at the position of the Sun with R = 8.5 kpc. The references of these 
curves are |106| in green, I107| in purple, |108| in blue and 1 1121 in red. 

2.4 Propagation models 

Various approaches are useful to solve the transport equation. The simplest ap- 
proximation is the so called "leaky-box" model (LBM) which was employed in 
pioneering studies. In the leaky-box approximation, the Galaxy is described as a 
finite and homogeneous volume with uniform gas density. Each nucleus escapes 
from this volume with a probability l/r esc . The diffusion term V ■ uDViVjJ can 

be expressed as —Ni/T esc . The LBM can be considered as a diffusion model in two 
limiting cases: 

• cosmic rays diffuse rapidly in the Galaxy and reflect at the Galaxy halo bound- 
ary with little leakage from the system; 

• the Galaxy halo is much flatter than the radius of the Galaxy with thin source 
and gas disks [114] . 

The leaky-box model has been used successful to explain most observed cosmic 
ray fluxes of stable nuclei, however, it cannot deal with the complexities such as 
spatially dependent source distributions, etc. More complete treatment of all the 
relevant processes and more realistic description of the Galactic ingredients use 
other techniques to solve equation |2 . 9| explicitly, in which two main ones have been 
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employed to date: analytical (or semi-analytical) models and purely numerical mod- 
els. Several software packages have been developed for this purpose. For example, 
USINE is a commonly used code employing the analytical approach but no public 
version has been released yet [1 151 I116| . GALPROP 117 is a publicly available 
numerical code used widely not only for cosmic ray nuclei, electrons but also for 
photons [93j [44] . An overview of these models will be presented in this section. 



2.4.1 Description of the Galaxy 

Any solution of the propagation equation is based on some fundamental assump- 
tions regarding the Galaxy, including its geometry, its matter content as well as the 
magnetic fields. Cosmic rays are thought to diffuse in some containment volume be- 
yond which they freely stream out. The density outside the boundary drops to zero. 
Radio observations of galactic halos indicate that the shape of the confinement vol- 
ume might radially follow the galactic disc, but with a greater thickness. Commonly 
the Galaxy halo is modeled with cylindrical symmetrically with radius R = 20 kpc 
and half-height z% whose value is still unknown but is reasonably believed to be 
greater than a few kpc. The density of cosmic rays satisfies the boundary condition 
N(r — R,z) — N(r, z — iz^,) = 0. The Galactic disk is embedded in the halo with 
half-height h ~ 100 pc. 



The gas density distribution 



As discussed in section [2. 1.2| and section 2.1.3[ cosmic rays may interact with inter- 



stellar gas causing secondary production of particles and energy losses. Therefore, 
the gas density is a basic ingredient to affect these processes. 

The ISM is a mixture of neutral atomic hydrogen HI, ionized hydrogen HII, 
molecular hydrogen H2 and helium components. The densities of these compo- 
nents vary with the radial distance r. In some cases, for example in the ana- 
lytical code USINE, it is taken as a simplified average gas density in the disk of 
~ 1 proton/cm 3 . In numerical code GALPROP, a more realistic gas distribution 
is used instead of a constant gas density. The H2 density is calculated from the 
CO volume emissivity [118] and the conversion factor from CO emissivity to nn 2 
is taken as 1.9 x 10 20 molecules cm -2 (K km s -1 ) -1 |112j . The HI distribution is 
taken from the model in |119j but renormalized to agree with the density distribu- 
tion perpendicular to the Galactic plane |120j and |121| . The HII distribution is 
taken from a cylindrically symmetric model [122J. The hydrogen number density 



distributions are plotted in figure 2.4 for height z =0, 0.1 and 0.2 kpc. The helium 



number density fraction in the gas is taken as 0.11 



The interstellar radiation field and the Galactic magnetic field 

The interstellar radiation field (ISRF) and the magnetic field have a strong influence 
on electron energy losses, 7-ray production from inverse Compton scattering and 
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Figure 2.4. The number density distribution for HI (dashed lines), HII (dotted 
lines) and H2 (solid lines) in the Galaxy, taken from 1441 . Curves for a specific 
gas type are arranged with decreasing density for z =0, 0.1 and 0.2 kpc (nu 2 at 
z = 0.2 kpc is too low to be shown in this figure). 



synchrotron radiation. These aspects are taken into account in numerical codes but 
can only be treated by assuming mean values in analytical codes. 

The Galactic interstellar radiation field (ISRF) results from emission from stars, 
and the scattering, absorption, re-emission of absorbed star light by dust in the 
ISM. Therefore, the estimation of the ISRF distribution is difficult and relies on 
the luminosity distribution from the stellar populations of the Galaxy, the dust 
distribution and the description of the absorption, scattering of star light and re- 
radiation processes. A recent calculation of ISRF can be found in |123l 124, 125 . 

The fine structure of the Galactic magnetic field is far from being fully under- 
stood. An assumption regarding the magnetic field is only used to calculate the 
electron synchrotron losses. A spatially dependent model adjusted to match the 
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408 MHz synchrotron longitude and latitude distributions j!26j . is used in GAL- 
PROP. 



2.4.2 Analytical approach 

Based on all the assumptions made on the Galaxy geometry, since the disk is much 
thinner than the halo size, the disk is considered as infinitely thin for practical 
purpose in the analytical approach. Cosmic ray sources and their interactions with 
the ISM are confined to this thin disk, and reacceleration is also assumed to take 
place in the disk. As a consequence, a factor 2hS (z) is added to the terms relating 
to the source and fragmentation processes, as well as the terms relating to energy 
losses and gains. The diffusion which occurs throughout the disk and the halo is 
assumed to have the same strength but does not have any spatial dependence. 
Assuming steady-state, the transport equation can be rewritten as a Laplace 



equation in a cylindrical geometry, by replacing V ■ 
tion|2~9|by 

" d 2 + 1 d_ 

dz 2 rdrir-. 
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The density can be obtained by solving the equation using a Bessel expansion 
method. One can expand all the quantities over the orthogonal set of Bessel func- 
tions [Jo (Cfcjj)] j n wn j cn k i s the order of the Bessel decomposition and Cfe 
are the successive zeros of function Jo: 

oo 

r 



k=l 



N{r,z)=Y,N k {z),h[C k -), (2.17) 



OO 

?(r)=53&J (Cfc^)- (2.18) 

k=l 

The solution of the cosmic ray density N (R, z) comprises contributions from 
the disk and from the halo. The contribution from the disk involves the primary 
sources and the spallation production concentrated in the disk as well as the energy 
losses and the diffusive reacceleration. The contribution from the halo involves the 
products from radioactive decay in the whole halo. For each contribution, the 
detailed expression of the solution can be found in [1151 1110] . 

The analytical approach has the advantage of showing a direct relationship be- 
tween propagation parameters. Another benefit is that the computation is fast. 
However, the analytical solution is only based on some simplified assumptions and 
thus cannot extend to more complicated cases. For example, anisotropic diffusion 
proposed in some literatures |127| is not able to obtain analytical solutions. It is a 
challenge to use analytical methods to treat electron energy losses and photon pro- 
duction since information on ISRF and magnetic fields are difficult to be imported 
in analytical codes. 
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2.4.3 Numerical approach 

The numerical solution of the transport equation is based on a implicit scheme 
in the GALPROP code. Terms such as diffusion, reacceleration, convection and 



energy loss in equation 2.9 can all be hnitc-diffcrcnccd for each coordinate (i?, z, 



p) or (x, y, z, p) in the form 

~dt " At At +qu (2 - 19) 

where all terms are functions of (R, z, p) or (x, y, z, p). 

To ensure stability for large time step At, the Crank-Nicolson method |128j is 
used in GALPROP, which is second-order accurate in time since in this method all 
the terms are alternatively finite-differenced in the form 

dN, N* +At - N! 



dt At 

2At + 2At 



(2.20) 



The detailed derivation and expression of coefficients a±, a-i and «3 for each term 
can be found in [53]. 

Compared to analytical programs, numerical approaches need a heavier compu- 
tation effort. Numerical methods have been developed to deal with both two and 
three dimensional spatial models, allowing us to handle more complicated and more 
realistic models involving spatially varying quantities, for example, anisotropic dif- 
fusion coefficients and electron energy losses. The production of photons can be 
treated more completely in numerical codes while only the hadronic production of 
photons can be dealt with in analytical programs until recently. 



2.5 Current status 

Investigations of cosmic ray propagation have been addressed using both analytical 
and numerical methods (e.g., most recently |1291 11301 WA 11311 [73]), applied to 
experimental data including stable and unstable nuclei, electrons and gamma rays. 
In a given propagation model, stable secondary-to-primary ratios can be used to 
determine the ratio of the halo size to the diffusion coefficient while the radioactive 
isotopes allow us to break the degeneracy between these two parameters. Moreover, 
the source spectrum can be accessed from the propagated fluxes of primary nuclei 
(mainly protons) and electrons. 

From the inference of cosmic ray isotropy and confinement, diffusion should 
inevitably be involved in the propagation processes. The existence of reacceleration 
and convection is not proved definitively. Therefore different models were studied 
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in literature, such as plain diffusion (PD) models, diffusion reacceleration (DR) 
models, diffusion convection (DC) models and diffusion reacceleration convection 
(DRC) models. No specific model can be considered as the best one so far to explain 
all the observations as well as being physically reasonable. The most favored model 



claimed in [1291 11301 l?2"] is the DRC model as shown in figure 2.5 It points to 5 
higher than 0.8 which is highly disfavored by the cosmic ray anisotropy problem, 
i.e. too large anisotropy is predicted compared to the observed one at highest 
energies > 10 14 eV. DR models can explain quite well the sharp peaks observed 
in the secondary-to-primary ratios (e.g., B/C, [Sc+Ti+V]/Fe) at energies around 
1 GeV/n but they cannot reproduce the proton and helium fluxes unless a break 
is introduced around 10 GeV in the injection spectra [44]. PD and DC models 
require a break in the rigidity dependency of the diffusion coefficient D, i.e. defining 
D as f3D (p/ po) Sl and as (3D (p/ p ) S2 below and above the reference rigidity po 



respectively, or/and an additional factor (3 V , as shown in figure 2.6 Both the break 



in the diffusion coefficient and the break in the injection spectra are arbitrary and 
not physically motivated. The factor f3 v , that only has an effect on non-relativistic 
particles, may be related to nonlinear MHD waves [103] . 




{GeV/n} 



Figure 2.5. B/C ratio compared with the ratio from best-fit DR (red) and DRC 
(blue) models. The shade areas are the 68% confidence level. Taken from 72 . 

Even in frameworks involving the same processes, the derived values of parame- 
ters vary in different studies. The published results do not always present consistent 
answers on the best-fit values of propagation parameters, especially the most rel- 
evant one, 8. For example, for PD models, the typical diffusion coefficient slope 5 
(means 62 if there is a break in S) obtained varies from 0.4 to 0.6 (e.g. [93" jll03|ll31j . 
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Figure 2.6. Left: B/C ratio for DC models with a break in the diffusion coefficient 
for convection velocity dV/dz = (solid lines), 5 (dotted lines) and 10 (dashed 
lines) km s~ 1 kpc -1 and z^=5 kpc. From |93| . Right: B/C ratio for a PD model 
with both a break in the diffusion coefficient and an additional /3 3 factor. Taken 
from 1X03] . 



For DR models, the best-fit S was determined to be about 0.3 in [53 0U [73] but 
about 0.5 in |131j and about 0.2 in |129[ll30j . The errors on S are usually not larger 
than ±0.05. 

Systematic uncertainties on quantities such as gas density, cross sections and 
data bias, are discussed in [130] , The gas density has a small impact on the diffu- 
sion slope but has very strong effects on other parameters like the normalization 
of the diffusion coefficient Do, the Alfven speed va and the convection velocity V c . 
Nuclear cross sections are crucially related to the destruction of primary cosmic 
rays and production of secondary cosmic rays. Their influences arc model depen- 
dent, however. The typical variance in the best-fit values is a factor of 2 for Dq, 
~10% for S, ^50% for va and ~5% for V c . Another important uncertainty arises 
from data bias since errors might be underestimated by a given experiment. The 
value of S determined by using data sets from different experiments can vary by 
more than 0.3 [1291 1130j . Therefore, to improve the reliability of parameter de- 
termination, more realistic gas density distributions and cross sections should be 
employed. Furthermore, the systematic uncertainties could also be reduced by us- 
ing spectra or spectrum ratios for all the necessary species provided by a single 
experiment to avoid data sets inconsistency arising from using data from different 
experiments. Therefore using PAMELA data exclusively to constrain propagation 
models is motived in my work and will be further discussed in chapter [5] 



Chapter 3 

The PAMELA experiment 



The PAMELA experiment is a satellite-borne apparatus designed to identify and 
measure charged particles and especially antiparticles in the cosmic radiation |63j . 
PAMELA is installed inside a pressurized container attached to a Russian Resurs- 
DK1 Earth-observation satellite which was launched into Earth orbit by a Soyuz-U 
rocket on June 15th 2006 from the Baikonur cosmodrome in Kazakhstan. The 
container of PAMELA is connected to the satellite body with a mechanical arm 
which can move the container from the parked position with downward orientation 
in which it is kept during launch to the position with upward orientation kept 
during data acquisition mode. 

Until now the instrument has been traveling around Earth along an elliptical 
and semi-polar orbit for almost six years, with an altitude varying between 350 
km and 600 km, at an inclination of 70 degrees. The trajectory thus goes through 
regions with varying geomagnetic cutoff, which effects the incident cosmic ray flux, 
and also passes the outer electron belt and the South Atlantic Anomaly (SAA) 
(figure 3.1 ). 

In this chapter, the scientific objectives of PAMELA will be presented in section 
13.11 and each sub-detectors of the instrument will be described in section 13.21 



3.1 Scientific objectives 

The design goal for PAMELA performance is to measure particle and nuclei fluxes 
over a wide energy range, and with unprecedented precision as a long exposure 
is achieved and no residual overburden of atmosphere needs to be compensated. 
Particularly, compared to previous experiments, PAMELA extends the energy range 
of antiprotons and positrons to both higher and lower energies. The statistics 
exceeds previous experiments by more than one order of magnitude after three 
years data taking. 

Table [3~T] shows the nominal design goals for PAMELA performance. 
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Rotational 




Figure 3.1. Orbit of the Resurs-DKl satellite. The satellite is travelling around 
Earth along an elliptical orbit, at an altitude ranging between 350 km and 610 km 
with an inclination of 70°. 



Energy range Statistics (3 years) 

Antiprotons 80 MeV - 190 GeV W 

Positrons 50 MeV - 270 GeV 10 5 

Positrons+Electrons Up to 2 TeV (from calorimeter only) 

Electrons 50 MeV - 400 GeV 10 6 

Protons 80 MeV - 700 GeV 10 8 

Light nuclei (up to Z=6) 100 MeV/n - 250 GeV/n He/Be/C 10 7 - 74 / 5 



Antinuclei search Sensitivity of O (10 7 ) for Antihelium/helium 

Table 3.1. PAMELA design goals. 

The PAMELA mission mainly focuses on the precise measurement of antiprotons 
and positrons. By studying the antimatter component of the cosmic radiation, the 
following themes will be addressed: 

• To search for evidence of dark matter particle annihilations by precisely mea- 
suring the antiparticle (antiproton and positron) energy spectrum. 

• To search for primordial antinuclei (e.g. antihelium) and to study low energy 
particles (e.g. trapped particles in the Earth's magnetic field and solar flare 
particles). 

• To test cosmic ray propagation models through precise measurements of the 
antiparticle energy spectrum and precision studies of light nuclei. 
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Besides, a reconstruction of the cosmic ray electron energy spectrum up to 2 TeV 
may give a hint for possible contribution from local sources. 

3.2 The detectors 

In order to reach the design goals of performance, PAMELA comprises several 
subdetectors, each providing an independent measurement of the incident particles. 
Figure |3.2| presents a schematic overview of the PAMELA instrument and shows 
the location of each subdetector. 



TOF(S1J |—p^ □□□□□□□□ < ]— ] 




Figure 3.2. A schematic overview of PAMELA instrument. The apparatus is 
~ 1.3 m tall, with a mass of 470 kg. Taken from |63j . 

The core detector of PAMELA is a 0.43 Tesla permanent magnet spectrometer 
(tracker) equipped with 6 planes of double-sided silicon detectors, allowing the sign, 



40 



Chapter 3. The PAMELA experiment 



absolute value of charge and momentum of traversing charged particles to be deter- 
mined. The spectrometer geometry and dimensions define the overall acceptance 
of the experiment which is 21.5 cm 2 sr. The maximum detectable rigidity is found 
to be ~ 1 TV from test beams. Spillover effects limit the upper detectable antipar- 
ticle momentum to ~ 190 GeV/c (~ 270 GeV/c) for antiprotons (positrons). The 
spectrometer is surrounded by a plastic scintillator veto shield which can be used to 
reject particles not cleanly entering the acceptance. An electromagnetic calorime- 
ter mounted below the spectrometer measures the energy of incident electrons and 
allows topological discrimination between electromagnetic and hadronic showers, 
or non-interacting particles. Planes of plastic scintillator mounted above and below 
the spectrometer form a time-of-flight system. It provides the primary experimen- 
tal trigger, identifies albedo particles, measures the absolute charge of traversing 
particles and also allows proton-electron separation below ~ 1 GeV/c. The volume 
between the upper two time-of-flight planes is bounded by an additional plastic 
scintillator anticoincidence system. A plastic scintillator system mounted beneath 
the calorimeter aids in the identification of high energy electrons and is followed 
by a neutron detection system for the discrimination of high energy electrons and 
hadrons which shower in the calorimeter but do not necessarily pass through the 
spectrometer. 

The PAMELA instrument is 1.3 m tall, with a mass of 470 kg. The average 
power consumption of PAMELA is 355 W, which is provided by the solar panels 
or batteries of the host satellite. Data are down-linked a few times per day to the 
mass memory of the satellite during acquisition, and radio-linked down to Earth 
when passing the ground center in Moscow, NTsOMZ. The average volume of data 
transmitted per day is about 15 GBytes, corresponding to ~ 2 million collected 
events. 



3.2.1 The Magnetic Spectrometer 

The magnetic spectrometer ,133! is designed to give a precise measurement of mo- 
mentum and charge (with sign) of the incident particle, as well as satisfying the 
requirements of the mission imposed by the satellite specifications. A compact 
mechanical assembly has been chosen and tested to withstand the stresses during 
the launch phase. The spectrometer is composed of a permanent magnet with an 
internal rectangular cavity, and a tracking system with six planes of double-sided 
silicon microstrip detectors, uniformly positioned along the cavity. Each plane in- 
dependently measures both the X and Y coordinates of the crossing point of an 
incoming ionizing particle. 

The reconstruction of the trajectory is based on the impact points and the 
resulting determination of the curvature due to the Lorentz force. The equation 
of motion describing a charged particle (with mass m and charge q) moving in a 
magnetic field B is: 
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d 2 f ( dr -\ . 

^ = q {d-t XB )> ^ 

where r is the position of the charge particle and 7 = 1/ \J\ — v 2 jc 2 . 

Introducing the path length I = /3ct and using p — m^j3c, equation (3.1 ) can be 
rewritten as: 

d 2 f q f dr ^\ f dr -\ , 

^ = ^r c {di xB ) = v{iu xB )> (3 - 2) 

where j3 = v/c and the magnetic deflection r\ is defined as the inverse of the 
rigidity R (R = p/q): 

„-I-f. (3,) 



Equation (3.2 1 can be solved by numerical methods, thus the deflection 77 of 
the particle is derived by looking for the set of initial conditions which best fit 
the measured coordinates of the particle trajectory. The spectrometer can also 
measure the absolute value of the charge since the ionization energy loss deposited 
in the sensitive areas of one plane is proportional to the square of the charge of the 
particle. 

The upper limit of the detectable energy range PAMELA can achieve is con- 
strained by the spectrometer bending power, expressed by the so-called Maximum 
Detectable Rigidity (MDR). MDR is defined as the measured rigidity which corre- 
sponds to 100% uncertainty. This feature conflicts with the detector acceptance, 
expressed as a Geometrical Factor (GF), which is defined as the factor of propor- 
tionality between the detector counting rate and the intensity to isotropic radiation. 
While the acceptance grows with the cross section of the cavity, the bending power 
improves for a longer cavity and a larger magnitude of B. Given a constraint on 
B, a longer magnetic cavity enhances the MDR while lowering GF, and conversely, 
a wider acceptance increases GF but worsens the MDR since it is more difficult 
to maintain a high field over a larger area. The geometric design should give the 
best compromise between these two features. Since extending the measurement of 
antiprotons and positrons to higher energy is a main objective, the MDR has been 
preferred in designing the spectrometer. The measurement error for momentum p 
depends on two contributions, the finite spatial resolution of the tracking system a 
and multiple Coulomb scattering of the particles crossing the spectrometer, whose 
relative weight varies with momentum. The errors from these two contributions, 
expressed as Ap res and Ap ms can be derived as: 

Apres o- 

— « m * (3-4) 

- 2 \ 2 

(3.5) 



and 
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where p is the momentum of the particle, L is the length of the track when pro- 
jected on the bending plane, and a is the spatial resolution. The measured spatial 
resolution with test beams is (3.0 ±0.1) /zm and (11.5 ±0.6) /xm in the bending 
and non-bending views, respectively [63 . Since the rigidity R = p/q, for a particle 
with a certain charge q, the relative error on rigidity is AR/R = A.p/p (see figure 



3.3 1. While spatial resolution plays a crucial role at high energy, multiple Coulomb 



scattering cause the main uncertainty of measured rigidity at low energy. Conse- 
quently, a long magnetic cavity with a strong magnetic field provides a good spatial 
resolution at high energy, while a minimal amount of material along the path of the 
particles reduce the scattering effect at low energy. The expected MDR is about 
1 TV/c and the computed GF for high-energy particles with straight tracks is about 
21.5 cm 2 sr [134]. 




Figure 3.3. Spectrometer resolution as a function of rigidity. The dotted lines 
present the rigidity relative error due to the finite spatial resolution AR res /R and 
the error due to the multiple scattering AR ms /R . The solid line shows the quadratic 
sum of the two. Taken from 135 . 

The upper rigidity limit for particles like protons, nuclei and electrons is directly 
connected to the MDR, but this is not the case for the antiparticles due to their 
rarity in cosmic rays. When the energy of particles increase, the tracks get closer 
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and closer to a straight line. The finite spatial resolution makes it difficult to 
properly determine the charge sign, which is used to distinguish antiparticles from 
particles. This effect, called spillover ; causes a non-negligible background when 
measuring antiparticles at high energy, especially for antiprotons since the number 
of protons is much larger than the number of antiprotons (about a factor of 10 4 at 
10 GeV). A simulation of proton spillover into a sample of antiprotons implies the 
detection of antiprotons is limited to about 190 GeV (see figure 3.4). 




10 ' l 10 u> 2 10 3 



Kinetic energy (GeV) 



Figure 3.4. Simulated spillover effect in the antiproton flux measurement. The 
points show the expected antiproton flux in three years measurements with PAMELA 
according to pure secondary production. The shaded area shows the simulated 
proton spillover in the antiproton sample. Taken from 1 1351 ■ 



The magnet 

The magnet is a 43.66 cm high tower, formed by five identical modules with a 
central rectangular cavity (16.14 cm x 13.14 cm). Each module is composed of 12 
Nd-Fe-B alloy elements with high residual magnetization (« 1.32 T). A picture of 
the entire configuration is shown in figure [375] 

This configuration has been chosen to have very high and uniform field strength 
inside the cavity and the lowest possible field intensity outside. To protect the 
magnetic material from chemical attacks, a 500 /mi aluminum layer covers all the 
free surfaces of the magnet. The field inside the cavity is almost uniform, with 
practically all the strength along the negative Y-direction [136] . As a consequence, 
particles are bent in the XZ plane within the cavity, due to the Lorentz force F = 
qv x B. The magnetic field has been mapped by means of an FW-Bcll Gaussmeter 
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equipped with a three-axis Hall probe mounted on an automatic positioning device. 
The measured main field component (B y ) is plotted for the central plane of the 
cavity and plotted for the central axis of the cavity in figure |3.6| In the center 
of the cavity (z=0) the value reaches 0.48 T and remains nearly constant across a 
wide region. The measured average magnetic field is about 0.43 T. 




Figure 3.5. The magnet tower (left). Sketch of a prototype of one magnet module 
(right). Taken from [T32] , 
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Figure 3.6. The main magnetic field component (B y ) plotted for the central plane 
(z = 0) of the cavity (left) and the (B y ) plotted as a function of z coordinate along 
the central axis (x = 0, y = 0) of the cavity (right). Taken from [132) . 

Silicon tracking system 

The tracking system is composed of 6 planes of high-precision silicon microstrip 
detectors, placed between the five magnetic modules and above and below the 
openings of the magnetic tower, with an uniform vertical spacing of 8.9 cm. Each 
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plane, housed in an aluminum frame, consists of 3 independent sections (ladders) 
along the X axis. Each ladder is formed by 2 rectangular double-sided n-type 
silicon sensors with dimensions 53.33 mm x 70.00 mm x 30 fim and a hybrid 
circuit which houses the front-end electronics (figure 3.7). On the junction side 
2035 p + type microstrips are implanted with a pitch of 25.5 ^m and on the ohmic 
side 1024 n + type microstrips are implanted with a pitch of 66.5 /xm. Strips on 
opposite sides are orthogonal and the spatial information of the impact point of the 
incident particle can be measured by looking at which strip collected the ionization 
charge on junction (X) view and ohmic (Y) view. 



Junction Side (X) 



8* 128= 1024 channels 
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Figure 3.7. The sketch of the strips layout on the junction side and ohmic side of 
the ladder. Both the read-out electrodes perpendicular to the n+ strips and their 
connections on the diagonal of the sensor are shown on the ohmic side. Taken from 

EH- 



3.2.2 The time of flight system 

The time of flight system (ToF) is designed to fulfill several goals |137j : 

• provide a fast signal for triggering data acquisition of the whole instrument; 
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• measure the flight time of particles crossing its planes; once this information 
is integrated with the measurement of the trajectory length through the in- 
strument, the particle velocity ft can be derived. This feature enable also 
the rejection of particles entering the apparatus from below, called albedo 
particles; 

• determine the absolute value of the charge Z of incident particles through the 
multiple measurement of the ionization energy loss dE/dx in the scintillator 
counters. 

Additionally, segmentation of each detector layer in strips can provide a rough 
tracking of particles, thus helping to reconstruct their trajectory outside the magnet 
volume. 



The ToF (see figure 3.8) is arranged in three planes, referred to as SI, S2 
and S3, and each composed of double layers of segmented plastic scintillators to 
improve the reconstruction efficiency for crossing particles. SI is placed on top 
of the experiment, with eight 330 mm x 51 mm paddles forming its first layer 
Sll and six 408 mm x 55 mm paddles form its second layer S12. The overall 
sensitive area of each layer of SI is 330 x 408 mm 2 . S2 and S3 are placed above 
and below the spectrometer respectively. The two layers of S2, S21 and S22, are 
both divided into two paddles, with dimension 180 mm x 75 mm for S21 and 150 
mm x 90 mm for S22, resulting in a sensitive area 150 x 180 mm 2 for each layer. 
For S3, the first layers S31 is divided into three 150 mm x 60 mm paddles, while 
the second layer S32 is divided into three 180 mm x 50 mm paddles. The overall 
sensitive area of each layer of S3 is 150 x 180 mm 2 . For each plane, the paddles 
of the upper layer are orthogonal to those of the lower layer, therefore allowing a 
two dimensional coordinate measurement of the impact points of charged particles. 
BC-404 manufactured by Bicron company was chosen for the scintillator material, 
characterized by a rise time 0.7 ns and a decay time 1.8 ns. The two ends of 
each paddle are read-out by a Hamamatsu R5900 PMT, which can achieve an 
amplification of about 4 x 10 6 at 900 V. Since the core of the PAMELA apparatus 
is a permanent magnet, all the PMTs have been shielded with a 1 mm thick /i-metal 
screen to avoid the influence of any residual magnetic field. 

The anode pulse of each PMT is converted both in charge and time, by connect- 
ing to an analog-to-digital converter (ADC) and a time-to-digital converter (TDC) 
respectively. When a charged particle crosses a layer, the ADC measures the ioniza- 
tion energy loss, and the TDC provides the relative time. The charge identification 
capabilities of ToF were evaluated during a test with particle beams performed at 
the GSI laboratory in Germany, which indicates the measured charge uncertainty is 
less than 0.1 for protons and 0.16 for carbon |138j . The combined TDC information 
of all the ToF planes is used to generate the main PAMELA trigger and determine 
the flight time of the incoming particle. The standard trigger configuration re- 
quires the coincidence of at least one TDC signal from each of the three planes. In 
the radiation belts and inside the SAA the requirement on SI is removed as SI is 
saturated by low energy particles. Moreover, TDC information can establish the 
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Figure 3.8. The ToF system. The sensitive areas are 330 X 408 mm 2 for SI and 
150 X 180 mm 2 for S2 and S3. The distance of SI and S3 planes is 77.3 cm. Taken 
from [132] , 



incident particle direction by checking the order in which the layers have been hit. 
This is important since the sign of charge can be determined by the deflection only 
if the direction of motion is known. 

The velocity j3 (v/c) of an incident particle can be calculated by measuring the 
time needed for the particle traveling from SI to S3. For a particle with velocity 
(3, momentum p and mass m, it follows that: 

P= , 1 =, (3.6) 
Jl+ (mc 2 /pc) 

which give the possibility to discriminate types of particles with different masses 
at low energy. The measured ToF time resolution about 250 ps allows electrons 
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(positrons) to be separated from antiprotons (protons) up to 1 GeV/c [63] 



3.2.3 Anticoincidence system 

The primary aim of the Anticoincidence (AC) system designed and built at KTH is 
to identify events yielding "false" triggers, which might be generated by secondary 
particles produced in the mechanical structure of the experiment. It also can help 
to reject out of acceptance events. 

The PAMELA experiment contains two AC systems. One of them consists 
of 4 plastic scintillators (CAS) covering the sides of the magnet each of which 
has an approximate rectangular shape, and 1 scintillator covering the top (CAT) 
which has a star shape with a rectangular hole in the center corresponding to the 
acceptance of the spectrometer (see figure 3.9). The other one consists of 4 plastic 
scintillators (CARD) surrounding the empty volume between SI and S2. All the 
scintillators are Bicron BC-448M and each CAS and CARD detector is read out 
by two identical Hamamatsu R5900U PMTs in order to decrease the possibility of 
single point failure, while the CAT detector is read out by 8 PMTs for the same 
reason and also to cover the irregular shaped area. The scintillators and PMTs 
are housed in aluminum containers which provide light-tightness, allow fixation to 
the PAMELA superstructure and ensure that a reliable scintillator-PMT coupling 
is maintained. No additional magnetic shielding is required due to the small fringe 
field from the magnetic spectrometer at the position of the PMTs. A particle 
traversing an AC detector is registered as a hit if it deposits at least ~ 0.8 MeV 
energy in the scintillator. The detection efficiency for charged particles is measured 
to be 99.9%. 




Figure 3.9. A schematic view of CAS (purple) and CAT (green). The CAS scintilla- 
tor is ~ 40 cm tall and 33 cm wide. The hole in the CAT scintillator is ~ 22 X 18 cm 2 . 
Taken from [T32j . 
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3.2.4 Electromagnetic Calorimeter 



The calorimeter is the key detector used to select antiprotons and positrons from 
like-charged backgrounds which are significant more abundant in cosmic rays. An- 
tiprotons must be separated from a background of electrons that decreases from 
~ 10 3 times the antiproton component at 1 GeV/c to less than 10 2 above 10 GeV/c, 
and positrons from a background of protons that increases from ~ 10 3 times the 
positron component at 1 GeV/c to ~ 5 x 10 3 at 10 GeV/c. This means that the 
PAMELA detectors have to separate electromagnetic from hadronic particles at a 
level of 10 5 — 10 6 . Most of this separation is provided by the calorimeter. 

The calorimeter is composed of 44 silicon sensor layers interleaved with 22 planes 



of tungsten absorbers (see figure 3.10). Each tungsten layer has a thickness of 26 
mm which corresponds to 0.74 X (radiation lengths), giving a total depth 16.3 
Xo or 0.6 A (interaction lengths). Each silicon plane consists of a 3 x 3 matrix 
of 8 x 8 cm 2 silicon detectors. The detectors are separated from each other by 
a ~ 35 mm wide non-sensitive area, called a dead area. About 5% area of one 
plane is covered by such dead areas. Each silicon detector is 380 fim thick and 
segmented into 32 strips with a pitch of 2.4 mm. Two consecutive sensor layers 
and one sandwiched tungsten absorber form a detector plane. The orientation of 
the strips of the two layers in a detector plane is orthogonal and therefore provides 
two-dimensional spatial information of a particle shower. 

The longitudinal and transverse segmentation of the calorimeter, combined with 
the measurement of the particle energy loss in each silicon strip, allows a high rejec- 
tion power of electrons in the antiproton sample and protons in the positron sample. 
A good agreement is found between simulated and experimental calorimeter data. 
Simulations demonstrate a rejection factor of about 10 5 for electrons in antiproton 
measurements with 90% antiproton identification efficiency |139| . The calorimeter 
is also used to reconstruct the energy of the electromagnetic showers, providing a 
measurement of the energy of the incident electrons independent from spectrome- 
ter, thus allowing a cross-calibration of two energy measurements. The calorimeter 
energy resolution has been measured as ~ 5.5% up to several hundred GeV (shown 



in figure 3.11). In order to measure very high energy electrons (~ 300 GeV to > 1 
TeV), calorimeter is equipped with a self-trigger capability. A self-trigger signal 
is generated when a specific energy distribution is detected predetermined planes 
within the lower half of the calorimeter. By requiring that self-triggering particles 
enter through one of the first four planes and cross at least 10 radiation lengths, the 
geometrical factor can achieve 600 cm 2 sr, which is about 30 times larger than the 
default PAMELA geometrical factor defined by the magnetic spectrometer. Since 
the geometrical factor is highly increased in self-trigger mode, PAMELA has the 
capbility to measure the very-high energy electrons which are rare in the cosmic ra- 
diation. The calorimeter energy resolution in self-trigger mode is estimated through 



simulation to be ~ 12% up to about 800 GeV, as shown in figure 3.11 
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Figure 3.10. The PAMELA electromagnetic calorimeter comprising 22 calorimeter 
modules. The device is ~ 20 cm tall and the active silicon layer is ~ 24 X 24 cm 2 in 
cross-section. The bottom plot shows the detail of a single module consisting of a 
tungsten layer sandwiched between two silicon detector planes. Taken from [63| . 



3.2.5 Bottom Scintillator S4 

The bottom scintillator S4, referred to as the shower tail catcher, is used to improve 
the PAMELA electron-hadron separation by measuring shower leakage from the 
calorimeter. S4, with a sensitive area of 482 x 482 mm 2 and a thickness of 10 
mm, is located directly beneath the calorimeter and read out by six PMTs placed 
along the two opposite sides. The S4 detector detects showers not contained in the 
calorimeter. When the signal in S4 exceeds 10 MIPs (where 1 MIP is the most 
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Figure 3.11. The calorimeter energy resolution as a function of the incident electron 
energy Ei n . The filled circles are for normal operation (experimental data) and the 
open circles are for the self-trigger mode (simulations). Taken from |63| . 

probable energy deposited by a normally incident minimum ionizing particle) and 
coincides with the main trigger signal, an on-board neutron detector is read out. 



3.2.6 Neutron Detector 

The primary purpose of the neutron detector (ND) is to complement the electron- 
proton discrimination capabilities of the calorimeter. Also, combined analysis of 
calorimeter and ND information will expand the energy range for detected primary 
electrons up to 10 TeV. 

The detector is located below the S4 scintillator and consists of 36 gas propor- 
tional counters stacked in two planes of 18 counters each, oriented along the y-axis 
of the instrument. The size of the neutron detector is 600 x 550 x 150 mm 3 . The 
counters arc filled with 3 He and surrounded by a polyethylene moderator enveloped 



i 



a. -a 

0.16 
0.14 ■ 

C.-2 
• 0.1 
0.08 
G.3B 
0.04 ■ 
0.02 

n 



52 



Chapter 3. The PAMELA experiment 



in a thin cadmium foil to prevent thermal neutrons entering the detector from the 



sides and from below (see figure 3.12) 




Figure 3.12. The neutron detector partially equipped with 3 He proportional coun- 
ters. The neutron detector covers an area of 60 X 55 cm 2 . Taken from |63| . 

When a high-energy hadron interacts inside the calorimeter, a large number of 
neutrons are produced by the decay of excited nuclei, while the number is 10-20 
times lower if the primary particle is an electron and the neutrons are generated 
via photo-nuclear interactions. A part of these neutrons is thermalized by the 
polyethylene moderator and detected by the 3 He counters. 



To summarize, this chapter details the scientific objectives of the PAMELA ex- 
periment and describes all the sub-detectors including the spectrometer, the time- 
of-flight system, the anticounters, the calorimeter, the bottom scintillator and the 
neutron detector. Using a combination of these detectors, different cosmic ray 
species can be identified and their fluxes measured. In chapter [4j cosmic ray an- 
tiprotons will be selected with the help of PAEMLA instrument and the flux of 
antiprotons in the cosmic radiation as well as the antiproton-to-proton flux ratio 
will be reconstructed . 



Chapter 4 



The antiproton flux and 
antiproton-to-proton flux 
ratio 

Cosmic rays are dominated by protons and helium nuclei, with a small fraction 
of electrons, helium-3, deuterium and heavier nuclei, as well as rare antiparticles. 
In the negatively charged part, antiprotons constitute only a small fraction (about 
10 -3 ) compared to the main component, electrons. Besides primary cosmic rays, 
PAMELA also records particles created in the interactions of primary cosmic rays 
with the experiment materials, for example positive and negative pions, which 
complicates the identification of antiprotons. 

The procedure used to determine the antiproton flux and antiproton-to-proton 
flux ratio (p/p) contains three steps: 

1. Select a reliable antiproton sample. Section [4?T] lists the selection criteria. 

2. Calculate the efficiencies of the selection cuts, i.e. the probability that an 
antiproton will pass the selection criteria. This is presented in section [472] 

3. Correct for other factors such as geometrical factor, hadronic interaction 
losses, the live time of measurements and transmission through the geomag- 
netic field. All these corrections are discussed in section l4~3l 

The final results of antiproton flux and p/p ratio measured by PAMELA are shown 
in section [L4l 

4.1 Antiproton selection 

Before antiproton identification, non-corrupted data (i.e. data contain proper in- 
formation from each sub-detector) are pre-selected. Events in the high radiation 
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region such as the South Atlantic Anomaly (SAA), yielding high counting rate 
which may cause unstable performance of detectors, are rejected. 

To obtain a clean antiproton sample, selection cuts have been applied to several 
variables. Charge one hadrons entering the instrument from above with a well 
reconstructed track in the tracking system and a trajectory contained in the fiducial 
geometric acceptance have been selected. Among the surviving particles, negatively 
charged particles are selected as antiprotons, while a MDR cut is applied to remove 
the oppositely charged contamination due to the spillover effect. Since the number 
of antiprotons is about 10 -4 of the number of protons, the MDR cut is crucial for 
the selection of rare antiprotons in order to reduce significant proton contamination 
in the antiproton sample. In following sections all the selection criteria are classified 
by sub-detector and are described in more detail. 

4.1.1 Tracker criteria 

The tracker cuts can be divided into four categories: (i) the basic tracker selection to 
reconstruct the track of the incident particle and to provide a reliable rigidity mea- 
surement; (ii) the geometrical selection to define the acceptance of the experiment; 
(iii) additional tracker cuts imposed to further clean the (anti) proton sample; (iv) 
a MDR cut to reduce the background of spillover protons in the antiproton sample 
at high energies due to the wrong assignment of charge sign. 

The basic tracker cuts 

The basic tracker selection cuts are: 

• A single physical track reconstructed by the track fitting algorithm. Unrea- 
sonable events for which the fit routine does not converge or x 2 (the goodness 
of the fit) is less than zero are excluded. Most multiparticle events will be 
rejected by this cut. 

• Number of hits on planes in the x-view and y-view: N x > 4, N y > 3, and the 
lever-arm (defined as the distance (in planes) between the upper and lower 
impact points) in the x-view > 4. This selection ensures a good quality of 
the track. The number of fit points is larger in the x-view than in the y-view 
since the rigidity construction is performed from the bending in the x-view. 
A larger number of fit points and a longer lever arm will give a better track 
reconstruction. 

• An upper limit on x 2 , which is a comparison between measured and recon- 
structed impact points in the tracker planes 

X 2 < 12.42 + 199.5 xt] 2 + 153.1 x ry 4 , (4.1) 



where 77 is the deflection defined in equation 3.3 Multiple tracks or particles 
suffering multiple scattering when they cross the tracker planes usually yield 
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rather high x 2 . The upper limit on the x 2 is chosen so that the efficiency of 
this cut is about 95% and constant with rigidity. 

The geometrical cut 

The geometrical selection requires that the particle track is contained inside the 
fiducial acceptance during its entire passage from SI to S3. A 1.5 mm margin is 
subtracted from each side of the 8 tracker planes and 6 ToF planes which defines the 
geometrical acceptance. Therefore, this fiducial acceptance is 92% of the nominal 



acceptance mentioned in section 3.2.1 This cut is necessary to exclude particles 
crossing the magnet walls which might give incorrectly reconstructed tracks, and 
to avoid efficiency underestimation caused by bad tracking. For high energies, the 
fiducial acceptance is about 19.90 cm 2 sr [140] . 

Additional tracker cuts 

Additional tracker cuts contain: 

• A cut on the ionization energy loss in the tracker to select singly-charged 
particles. The mean values of the dE/dx measurements on the 12 tracker 
planes are contained inside the region bounded by the red lines shown in 
figure |4~Tj This cut helps to exclude multiply charged particles, multi-particle 
events and particles interacting in the tracker material. 

• A cut to further remove multiple tracks in the tracker (hereafter referred as 
TrkMultipleTracksCut). Events which fulfill one of the following conditions 
are excluded: 

(a) apart from the track, at least 3 hits located on the same side of the track 
in the x-view and 3 hits totally in the y-view; 

(b) apart from the track, at least 2 hits located on the same side of the track 
in the x-view and 2 hits totally in the y-view when hit PMTs are not 
associated to the track passing through SI and S2. Information from the 
ToF is used since if there is at least one hit PMT outside the track on SI 
and S2, this indicates a higher possibility of inelastic reaction occurred 
above the tracker and a lower limit should be put on the number of hits 
apart from the track. 

• A x 2 cut imposed only for low energy particles: 

X 2 < 5.99 + 131 x t] 2 + 99.09 x ry 4 (4.2) 

This selection, which is about 90% efficient, places a stronger limit on the x 2 
than the one used in the basic tracker cuts. As multiple scattering effects are 
significant at low energy, this cut is conservatively applied below 14.6 GV to 
reject particles scattering in the tracker system which might cause an unreli- 
able track reconstruction. 
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Figure 4.1. Tracker dE/dx selection. Particles within the red lines are selected as 
proton (antiproton) candidates. 



• Requirements on high energy particles to clean the tracking position mea- 
surements, including a x 2 limitation, a cut placed on the maximum energy 
release on the first tracker plane and the maximum multiplicity to remove 
events with accompanying hits due to delta ray emissions, no bad strips and 
a spacial resolution less than 0.01 mm. 



The MDR cut 

A further MDR cut is imposed due to the finite spectrometer resolution. The 
MDR, which is evaluated for each event during the fitting procedure, is required 
to be larger than Cx the upper limit of the rigidity bin containing the rigidity of 
the event, where C represents a coefficient. Since MDR = I/A77, this cut allows 
to reject events with large associated deflection errors and significantly eliminates 



the spillover protons, as shown in figure 4.2 A coefficient of 10 is enough to 
remove all the spillover protons [141 j . However, since the estimated MDR is about 
1 TeV, C = 10 restricts reconstructed energies to less than 100 GeV. Therefore, 
a coefficient of 6 is chosen to compromise the rejection power of spillover protons 
and the upper limit of detectable energy. Moreover, "sub-bins" are introduced in 
highest bins to increase the statistics for high energy antiprotons , i.e. each of the 
last three rigidity bins is divided into 20 sub-bins. A cut, MDR > 6x the upper 
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limit of the sub-bin containing the rigidity of the event, is used instead of original 
one in the last three bins. In the final calculation, the spillover contamination has 
to be precisely determined. 




Figure 4.2. The MDR distribution as a function of deflection. The red line is 
MDR = 6 X rigidity and the blue line is MDR = 10 X rigidity. The spillover protons 
can be clearly observed as the dense area which spills into the negatively-charged 
side. Events below the curve are rejected. 



4.1.2 ToF system criteria 

The ToF system measures the time-of-flight of the incident particles thereby provid- 
ing a velocity determination. This information can be used to reject albedo particles 
which will cause an upward going proton to be misidentified as a downward going 
antiproton. It also helps to identify low energy (anti)protons. The dE/dx mea- 
surements can exclude heavier particles, low energy electrons and pions. Combined 
with the hit information in the two top ToF scintillators, multiparticle events or 
interactions above the tracker can be rejected. The detailed selection cuts are as 
follows: 

• Events satisfying following requirements are selected to remove multiparticle 
events where particles traverse different paddles on the same scintillator: 

(a) no more than 1 hit paddle on Sll, S12, S22, S21; 

(b) at least 1 hit paddle on SI and S2; 
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(c) no more than 2 hit PMTs outside the reconstructed track on Sll, S12; 

(d) if there is a hit paddle on Sll, S12, S21 or S22, its PMTs must be 
associated to the track extrapolated from tracker or there must be TDC 
signals belonging to that hit paddle. 

No cut is applied on S3 as the particles interacting below the tracker whose 
rigidities have already been proper measured should be selected in the sample. 

The y measurement given by the TDC signals on Sll or S22 must be within 
a 6 cm tolerance margin around the y coordinate extrapolated to that plane 
from the reconstructed track. This cut checks the consistency between the 
track measured by the TOF and the one measured by the tracker. 

Events with /3 consistent with the expectation within 5a for (anti)protons are 
selected to reject other kind particles up to a few GeV, as shown in figure [473) 
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Figure 4.3. distribution as a function of rigidity. Particles within the red lines 
are selected as (anti)protons candidates. 

dE/dx cuts in SI and S2, similar to those used in the tracker, are applied to 
select singly-charged particles. 

Events with no more than 1 hit PMT outside the track on Sll and S12 and 
no large energy release in Sll, S12, S21 and S22 are selected at low rigidities 
(below 14.6 GV) to reject particles interacting before the tracker, which are 
mainly pion background. 
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4.1.3 Anticoincidence criteria 

The anticoincidence system (AC) is very useful to reduce events interacting inside 
the apparatus and producing secondaries. As mentioned in section |3.2.3[ these 
kind of events are always accompanied by multiple particles and can yield "false" 
triggers. In order to remove them, requirements are placed on AC, as shown below: 

• Events with no signal in CARD scintillators arc selected; 

• Events with no signal in CAT scintillator are selected. 

No cut is put on the CAS scintillators since particles backscattered from the ca- 
lorimeter can potentially be registered by the CAS scintillators but should not be 
rejected |142j . 

4.1.4 Calorimeter criteria 

The main task of the calorimeter is to identify antiprotons from an electron back- 
ground which is significantly more abundant. As discussed in section ^. 2. 4[ the lon- 
gitudinal and transverse segmentation of the calorimeter, combined with a dE/dx 
measurement in each silicon strip, allows the rejection of electromagnetic showers. 

For an electron traversing matter, the radiation loss (so-called Bremsstrahlung) 
exceeds the collision loss above a few tens of MeV and dominates the energy loss 
as the energy increases. A simplified description of an electromagnetic shower as- 
sumes an incident electron of energy Eq will lose half its energy to a bremsstrahlung 
photon after one radiation length, Xq, which corresponds to about 2 planes of the 
calorimeter. Each bremsstrahlung photon will after one radiation length produce 
an electron-positron pair, which in turn radiates another photon. This multiplica- 
tion proceeds to a maximum depth until the energy of the produced secondaries 
reaches the critical energy E c below which ionization losses start to dominate. The 
maximum depth is given by t max = ln(Eo/E c )/hi2 radiation lengths. Thus in the 
cascade the number of particles rises exponentially to a broad maximum and after 
that the shower decays slowly. The lateral spread of the shower depends mainly 
on the longitudinal depth and does not significantly depend on the energy of the 
primary electron. Multiple scatterings in the absorber have an important effect 
on the lateral spread, yielding two components in the shower: a narrow, strongly 
collimated central part due to the high-energy particles depositing most of the in- 
cident energy and a peripheral component spreading out as the shower penetrates 
deeper and low energy particles are created. The transverse spread is measured 
in a unit called the Moliere radius (Rm), defined as the average lateral spread of 
an electromagnetic shower initiated by an electron of energy E c when the electron 
traverses one Xq of material. The electromagnetic shower is about 95% laterally 
contained in 2Rm, which is about 1.8 cm (7.5 strips) for tungsten. 

Unlike electromagnetic showers, a hadronic shower results from different in- 
elastic hadronic interactions and consists of a wide variety of particles such as 
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pions, protons and neutrons, with large fluctuations in multiplicity and energy loss 
between individual showers. On average, about half of the incident energy is ac- 
quired by the particles produced in inelastic interactions. The resulting secondaries 
thus have large transverse momentum and the hadronic shower tends to be more 
spread out laterally than an electromagnetic one. The longitudinal development of 
a hadronic cascade is described in units of (nuclear) interaction length, Xa, which 
has a value of 9.6 cm for tungsten at high energies, indicating that most incident 
hadrons will interact deeper in the calorimeter or even traverse the calorimeter 
without interacting. 

Examples of an electromagnetic shower and a hadronic shower are illustrated in 



figure [4~4| In order to separate hadrons and leptons, an energy dependent criteria 
has been developed based on both simulations and tests with particle beams 139 j. 
Several calorimeter variables are used for the lepton/hadron separation. The total 
energy deposit in the calorimeter, referred to as qtot, allows a powerful separation 
between leptons and hadrons. For electrons, since most showers are contained in 
the calorimeter, q to t is usually normally distributed for a given incident energy. For 
hadrons, the distribution of q to t is flat with a sharp peak at low energies for non- 
interacting particles. Therefore an energy-momentum match exists for electrons, 
i.e. the qtot /rigidity satisfies a quasilinear relation, while for hadrons qtot /rigidity 
assumes a lower value. Other calorimeter variables used in the analysis are described 
below. 



The starting point of the shower 

While a hadronic shower has a roughly uniform probability to start in any plane of 
the calorimeter, an electromagnetic shower is more likely to start in the first three 
planes. A variable used to characterize this difference is referred as noint, given by 

2 22 

noint — 2_, % ' h (4-3) 

j=i i=i 

where 9+j = 1 if the ith plane of the jth view has strips registering energies compat- 
ible with a minimum ionizing particle within 4 mm from the reconstructed shower 
axis, otherwise 8^ = 0. The variable noint will increase as the interaction starts 
in deeper planes. Therefore it assumes low values for electromagnetic showers, and 
takes higher values for a non- or partially-interacting hadron. The distribution of 



noint from flight data and simulation are shown in figure 4.5 

Another variable which is sensitive to the starting point of the shower is the 
ratio of energy deposited in a cylinder of diameter 2 strips {q pre sh) and the number 
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of strips hit in the same cylinder (np res h), i.e. q P resh/n pres h, shown in figure 
As an electron interacts immediately in the first planes, the average energy deposit 
in each strip is expected to be lower compared to a hadron which starts to interact 
deeper in the calorimeter or does not interact. 
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Figure 4.4. An event display of a 50 GeV electron (left) and proton (right) recorded 
at CERN SpS facility (taken from 1 1391 ^ ■ In the top part of the figure, the two views 
(X and Y) of the six silicon planes are shown inside the magnetic cavity. In the 
bottom part of the figure, the two views (X and Y) of the calorimeter are shown. The 
color scale shows the detected energy in each strip. Orange (blue) area corresponds 
to a high (low) energy deposit. Vacancies in some layers are due to that some strips 
were not used in the test. Evident topological and energetic differences between 
electromagnetic and hadronic showers can be seen. 
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Figure 4.5. The distribution of the variable noint. Top: noint derived from flight 
data, both for negatively charged particles (mainly electrons) and positively charged 
particles (mainly protons). Bottom: noint derived from simulated electrons and 
antiprotons. The events above the red line are selected as (anti)protons. 
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Figure 4.6. The distribution of the variable qpresh/npresh derived from flight 
data. The negatively charged particles are mainly electrons, and the positively 
charged particles are mainly protons. The events above the red line are selected as 
(anti)protons. 



The longitudinal and topological profile 

While the energy deposit of an electromagnetic shower decreases after the shower 
maximum and spreads out laterally, the hadronic showers deposit their energy ap- 
proximately uniformly and any maximum lies deeper in the calorimeter. A quantity 
related to the longitudinal profile is defined as 

2 plmax 

q CO re = ^2 ^2 Qhittj ■ i, (4.4) 
j=i i=i 

where Qhit^ is the energy released in the jth view of ith plane within a cylinder of 
radius 2Rm centered on the shower axis, and pl m ax is the calculated electromagnetic 
shower maximum for a given incident energy provided by the tracking system. 

Furthermore, for an electromagnetic shower, before achieving the shower maxi- 
mum, the shower multiplication is expected to increase with shower depth and the 
shower particles should be collimated along the shower axis. A variable called n core 
is used to reveal this behavior, given by 

2 plmax 

n C ore = ^2 ^ NhiUj ■ i, (4.5) 

3=1 i=l 
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where Nhitij is the number of hit strips in the jth view of ith plane within a 
cylinder of radius 2Rm centered on the shower axis, and the plane number pl max 
is closest to the calculated electromagnetic shower maximum of the jth view. 

The energy density in the shower core weighted by the depth in the calorimeter, 
Qcore/n core , which is more sensitive to the shower difference between hadrons and 



electrons, is finally used to separate hadrons and leptons as shown in figure 4.7 



The lateral profile 

While the lateral development of a hadronic shower depends on the traverse mo- 
mentum of produced secondaries, which usually carry about half of the incident 
energy and thus cause a wide lateral spread, an electromagnetic shower has a lat- 
eral spread due to the multiple scattering of low energy particles which is usually 
less broad. A variable called n cy i, which is the number of strips hit in the cylinder 
of radius 2Rm around the shower axis, is sensitive to the lateral profile. Due to the 
difference between the behaviour of leptons and hadrons in the calorimeter, n cy i 
assumes a higher value for electrons than for hadrons (shown in figure 4.7 as well). 



4.1.5 Selecting Galactic particles 

As discussed in section fl.l| the Earth's geomagnetic field prevents low energy par- 
ticles from reaching the atomosphere. In order to select Galactic cosmic rays, only 
events with rigidities larger than the minimum value needed for a cosmic ray to 
penetrate the geomagnetic field and reach PAMELA are selected: 

rigbin lowerlimit > cutott PA MELA = 1.3 x SVC, (4.6) 

where rigbin lowerlimit means the lower limit of the rigidity bin. The Stoermer 
vertical cutoff (SVC), as defined in section is estimated using the satellite 
position. A coefficient of 1.3 is used here to ensure a robust selection of Galactic 
particles. 
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Figure 4.7. The distribution of the variable q C ore/n core (top) and the variable n cy i 
(below). The negatively charged particles are mainly electrons, and the positively 
charged particles are mainly protons. The events below the red line are selected as 
(anti)protons. 
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4.2 Selection efficiencies 

In order to determine the antiproton flux, the number of antiproton candidates 
surviving all the selection criteria is compensated for selection efficiencies. The 
efficiency of a set of selection cuts can be calculated as the surviving fraction of 
events when applying those cuts to a selected sample of antiprotons. Efficiency 
samples can be obtained in different ways: using simulations, test beams or flight 
data. The purity of the sample can be controlled well by using the first two methods. 
However, conditions can change during the flight which make the test beam data 
less useful. For example, during PAMELA flight, the performance of the Viking 



VAl chips (see figure 3.7), responsible for the readout of tracker, has degraded 
with time. In some periods, a Sll PMT was not operational and the PMT high 
voltage levels of the ToF varied. These conditions are impossible to reproduce in a 
test with particle beams. They also complicate the simulation model. Therefore, 
the flight data itself is used which intrinsically include the detector performance 
over time. Since independent detectors have to be used to select the sample of 
particles and determine their rigidities, biased samples might be introduced if the 
response of different detectors are correlated. To better understand and estimate 
the efficiencies, simulations are therefore used to cross-check the results. 

As it is impossible to select an unbiased, statistically significant antiproton 
sample from flight data, the efficiencies are derived from a proton sample with the 
assumption that protons and antiprotons behave identically in the detectors except 
for inelastic interactions. Therefore the only exception is the calorimeter, where 
the antiproton efficiency is corrected from simulations, taking into account that 
antiprotons and protons have different cross section for inelastic interactions. The 
efficiencies of the selection criteria discussed in section [4. 1| are presented below, in 
the order they were applied. 



4.2.1 Basic tracker cuts efficiency 

The basic tracker cuts are used to reconstruct a good quality track and give a reliable 
rigidity. Thus the efficiency sample of basic tracker cuts should be obtained without 
using the tracker system. In the low energy region, the rigidity can be obtained 
with the ToF system as /3m/ \fl — /3 2 . Due to the limited time resolution of the 
ToF system, this method is only applied below 1.7 GV/c [143j . For higher energy 
particles, the ToF rigidity resolution worsens since a small difference in the time 
of flight measurements produces a large relative difference in the reconstructed j3. 
Therefore, the tracker system is the only detector which can determine the rigidity 
of particles for the rigidity range of interest in this work. Fortunately, the efficiency, 
depending on the energy released in silicon planes, the curvature of the track and the 
multiple scattering, is expected to be constant for relativistic particles. For particles 
with a certain charge, the energy deposited in silicon planes is proportional to /3~ 2 . 
A particle deposits energy in the silicon planes and consequently creates a number 
of clusters, where a cluster is defined as one or more strips in the sensitive silicon 
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plane with a signal 7 standard deviations from the intrinsic noise of the channel. 
The number of created clusters naturally depends on the amount of energy released 
in the plane. For low velocity particles, as the multiplicity of clusters increases while 
the velocity decreases, the probability that the tracking algorithm will find a unique 
track decreases. However, for relativistic particles where (3 approaches unity, the 
probability is expected to be constant and thus the tracker efficiency also. The 
efficiency effected by the curvature of the track is also constant for relativistic 
particles since they produce nearly straight tracks. The last effect, characterized as 
the scattering angle (6 ~ z/ (p/3)) where z represents the particle charge and p the 
momentum, is negligible for singly-charged and f3 ~ 1 particles. Hence, for the basic 
tracker efficiency calculation, a sample is selected without rigidity determination. 
Once the basic tracker efficiency is obtained, the basic tracker cuts can be used 
to select efficiency samples for other selection criteria, allowing rigidity dependent 
efficiencies to be estimated. 

Selecting an experimental proton sample 

To select a clean proton sample, a set of cuts are applied to raw flight data as 
follows. In particular, in order to ensure the tracks of incident particles are inside 
the fiducial acceptance without using the spectrometer, the calorimeter is used to 
identify particle trajectories. 

Single particle selection The anticoincidence criteria described in section [4. 1.3| 
is applied. The number of hit paddles on Sll, S12, S21, S22 is required to be 
not larger than one. The number of hit paddles on both SI and S2 should be 
at least one. 

Charge one particle selection The energy released on SI should be less than 
1.8 MIP. 

Downgoing and high energy particle selection 



• 0.92 < f3 < 1.10. This cut selects downgoing particles whose velocity 
is positive and relativistic particles with (3 close to 1, considering the 
resolution of (3 is about 0.08. 

• Particles must cross both the last x-view plane and the last y-view plane 
of the calorimeter. 

Geometry constraint Particle tracks are reconstructed inside the calorimeter it- 
eratively. At each step, the hits in the calorimeter are fitted to define a single 
track, and the most distant points are rejected. The same fitting procedure 
is repeated until no hit can be found departing from the track further than 
a certain distance. The extrapolated tracks with good x 2 are required to be 
inside the acceptance, defined by the geometry of the 8 tracker planes and 
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the 6 ToF planes with a 0.7 cm tolerance applied to each side of every plane. 
This is the minimum tolerance which guarantees no underestimation of the 
efficiency due to uncertainties in of the calorimeter track fit and mechanical 
tolerances [140] . This cut also reject pions which are low rigidity particles 
with highly curved trajectories and then enter the calorimeter with an in- 
clined track. Their tracks back-propagated from the calorimeter therefore are 
not expected to be inside the acceptance. 

Electron rejection The energy deposited in the calorimeter strip closest to the 
track divided by the total energy deposit in calorimeter should be larger 
than 0.8 (this cut is referred to as CaloNotlntCut). This cut discards all 
the particles interacting in the calorimeter which naturally includes electrons 
since they interact immediately when they traverse the calorimeter. 



The efficiency of basic tracker seiection 

After applying the cuts described above to flight data, a proton sample is selected. 
The efficiency is defined as the fraction of events in the sample passing the basic 
tracker cuts detailed in section [4. The efficiency calculated for each day on orbit 
is plotted in figure |4.8[ The observed decrease is due to the tracker degradation, 
i.e. the efficiency falls gradually as the number of malfunctioning VA1 chips in- 
creases. Some other tracker performance issues can cause a short-term change of 
efficiency. For example, in 2006, from October 7th to October 11th (corresponding 



day numbers are 89 to 93 in figure 4.8), the power supply for the tracking system 
had a problem hence the tracking system was working only for a very short time. 
As a result, efficiencies for those days are zero or close to zero. In order to solve 
this problem, the power supply was changed to a redundant backup and DSlQ 7 
was switched off until October 23rd, which means one layer was not working on the 
y view with a resulting drop in efficiency. The days where the tracker was switched 
off will be excluded from the flux reconstruction to avoid an overestimation of the 
live time and an underestimation of efficiency. 

The simulated basic tracker efficiency for the tracker configuration with mal- 
functioning VA1 chips as in flight during July 2006 is 91.3% with a variation less 



than 0.7%, as presented in figure 4.9 The simulated result is slightly higher than 
the value obtain from flight data for that month, which is (90.6 ±0.1)%. The 
discrepancy between simulation and flight data is expected because the simulation 
can not give complete information regarding the \ 2 - However, the result shows a 
rather constant efficiency (variation less than 0.7%), which is consistent with the 
expectation that the efficiency is independent of rigidity. 

The efficiency derived by selecting protons not interacting in the calorimeter 
might be biased since the interacting protons may cause particles to be back- 
scattered from the calorimeter and reduce the probability for the tracker algorithm 



1 Digital Signal Processor - a part of the tracker electronics. 
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Figure 4.8. Time evolution of the basic tracker cuts efficiency. 
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Figure 4.9. The simulated basic tracker cuts efficiency for July 2006. 
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to find the correct single track. However, this underestimation is estimated to be 
negligible [144] . 



4.2.2 Additional tracker cuts efficiency 

Since the efficiency of the basic tracker cuts has been estimated, these cuts can 
be used to select efficiency samples for the other selection cuts, the efficiencies 
of which can be referred to as relative efficiencies. For example, if the selection 
criteria applied to select candidates of a certain species are grouped as A, B, C, 
etc, the sample for C's efficiency should be obtained by applying A plus B and 
other necessary cuts to reject extra background. Here C's efficiency is a relative 
efficiency, while A's efficiency is called an absolute efficiency. 

Instead of checking that tracks extrapolated from the calorimeter fall within 
the acceptance, the efficiency sample for the additional tracker cuts is derived by 
considering tracks reconstructed by the tracker system, i.e. applying the basic 
tracker cuts and the geometrical selection described in section [4. 

The rigidity dependent efficiency of the additional tracker cuts, integrated over 
the whole live time, is shown in figure 4.10| (black points). A discontinuity occurs at 



14.6 GV since different cuts are applied below and above this rigidity. However, for 
the rigidity range where the same cuts are employed, the dependence on rigidity 
is fairly constant. The slight decrease at high rigidities may be caused by delta 
rays. As the energy increases, more delta rays are generated thereby releasing 
more energy in the silicon planes and producing a higher hit multiplicity. These 



events are subsequently rejected by the delta ray cut (explained in 4.1.1 ). The time 



variation of the overall efficiency is shown in figure |4TTT) The changes of the tracker 
performance are visible in the efficiency. 

A correlation between the ToF and the tracker exists, as the cut used to remove 
multiple tracks in the tracker, named TrkMultiplcTracksCut, is based on the PMT 
information from SI and S2. Therefore a cut using information from S2, which 
might be sensitive to TrkMultipleTracksCut, is removed from the sample selection 
criteria to derive the efficiency. The resulting efficiency is compared with the one 



produced using the S2 cut in figure 4.10 Evidently, the difference between the two 



cases is negligible and will be omitted in the calculation of the efficiency. 



4.2.3 ToF efficiency 

A sample of protons has been derived without using the ToF system. Single down- 
going charge one protons are selected with good quality tracks inside the fiducial 
acceptance. The basic tracker cuts and the additional tracker cuts are both used 
to choose a sample of singly-charged, single particle events. The tracks of these 
particles must be inside the fiducial acceptance. The AC cuts are used to further 
reduce secondary events. Since the deflection of particles is required to be posi- 
tive, the remaining background are electrons, positrons and albedo singly-charged 
particles. Therefore the cut CaloNotlntCut is applied to reject the electrons and 
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Figure 4.10. Additional tracker efficiency integrated over the whole live time using 
flight data. The black (red) points are derived by using (without using) cuts on S2. 
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Figure 4.11. Time evolution of the additional tracker cuts efficiency. 
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positrons. Events with positive deflections are required which means the remain- 
ing contamination is upward-going antiprotons and therefore negligible (~ 10~ 4 of 
protons). 
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Figure 4.12. ToF efficiency based on flight data and integrated over the whole live 
time. 



The rigidity dependent efficiency of the ToF cuts, integrated over the whole live 
time, is shown in figure [4T2| A discontinuity occurs at 14.6 GV since different cuts 
are applied below and above this rigidity. However, for the rigidity range where 
the same cuts are employed, the efficiency does not vary more than 1%. The time 
variation of overall efficiency is shown in figure 4.13 A fluctuation can be observed 
which is caused by a known variation in the ToF performance. For example, around 
the day number of 20, 40, 180, drops of efficiency correspond to changes in the PMT 
high voltage settings. Around the day 70, a failure of one Sll PMT introduces a 
decrease of efficiency. A variation around the 200th day is due to a change in TDC 
threshold. 

The events which produce delta rays above the spectrometer, or are backscat- 
tered from the calorimeter while not interacting before S3, should be part of the 
proton sample. However, a "no hit" requirement on CAT and CARD removes a 
fraction of this kind of events, resulting in an overestimation of the ToF efficiency. 
In order to estimate the correlation between ToF and AC, the simulated efficiencies 
are derived by using the AC cuts and without using the AC cuts, called Stof_ac and 



ttofjnoac respectively. The results are shown in figure 4.14 As expected, et f_ac 
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day number 



Figure 4.13. Time evolution of the ToF efficiency. 



is higher than £tof_noac, since AC cuts remove a fraction of particles which should 
be included in the sample and will be rejected by the ToF cuts. The difference 
between the two efficiencies is less than 0.5%. This is added to the ToF efficiency 
as a systematic error. 



4.2.4 Anticoincidence efficiency 

After applying the basic tracker cuts, the additional tracker cuts and the ToF cuts 
to flight data, the proton sample for estimation of the AC cuts efficiency is further 
cleaned by using the calorimeter cuts described in section |4.1.4| The AC efficiency 
is the fraction of events yielding a signal on CARD or CAT in the proton sample. 
Since all the particles interacting hadronically before S3 are removed in the sample, 
the particles hitting CARD or CAT are either delta rays produced above the tracker 
or particles back-scattered from the calorimeter. These kinds of particles are good 
candidates but are rejected by the anticoincidence system, therefore the selected 
candidates need to be corrected for this effect. The derived efficiency is about 96% 
over the entire rigidity range, decreasing as the rigidity increases, as shown in figure 
|4.f 5| There is no significant time dependence since the AC performance is stable 
(see figure 4.16). However, a small drop can be observed correlated to the ToF 
performance. This drop is caused by the ToF system which has a lower rejection 
efficiency for delta rays or backscattered particles when selecting a proton sample. 
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Figure 4.14. Simulated ToF efficiency. The red (black) one is the efficiency derived 
by using (without using) AC cuts. 



Thus more events in the sample will give a signal in the anticoincidence system, 
yielding a lower AC efficiency during those days. The correlation between ToF and 
AC has already been discussed in section [4. 2. 3[ 



4.2.5 Calorimeter efficiency 

The calorimeter efficiency is derived for flight protons by calculating the fraction 
of events passing the calorimeter criteria over the sample obtained by applying 
the basic tracker cuts, the additional tracker cuts, the ToF cuts and the AC cuts. 
Since only positive deflection events are selected, the sample should only consist 
of positively charged particles, i.e. protons and a negligible number of positrons 
(~ 1CT 3 of protons). 

Antiprotons and protons behave differently in the calorimeter mainly due to 
their interaction cross sections. At 2 GeV, the cross section for a pp interaction is 
about 40 mb larger than a pp interaction |145] , The difference decreases at high 
energy but can not be assumed approximately equal until approximately 100 GeV. 
Therefore a systematically lower calorimeter efficiency for antiprotons than for pro- 
tons is expected until ~ 100 GeV. This is consistent with the simulated antiproton 
and proton efficiencies as presented in figure |4~T7| Hence, the antiproton efficiency 
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Figure 4.16. Time evolution of the AC efficiency. 
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is estimated by scaling the flight proton efficiency by the ratio of simulated antipro- 
ton efficiency and proton efficiency. The scaling factor is derived from simulation, 
however, this factor may not be exactly equal to the discrepancy appearing during 
flight. To account for the discrepancy between flight efficiency and simulated an- 
tiproton efficiency, a conservative systematic uncertainty of 5% is assigned to the 
antiproton calorimeter efficiency in the final calculation. 




Figure 4.17. The calorimeter efficiency for protons evaluated from flight data 
(black) and simulation (green). The blue points show the calorimeter efficiency for 
antiprotons from simulation. The red points are the calorimeter antiproton efficiency 
scaling the flight proton efficiency by the ratio of simulated antiproton efficiency and 
simulated proton efficiency. 



4.2.6 MDR efficiency 

The MDR efficiency is derived from the sample obtained by applying all the other 
cuts implemented in section [4 . 1 1 apart from the MDR cut. As mentioned in section 
|4.1.1| in the last three bins the MDR cut requires a MDR larger than 6 times the 
upper limit of the sub-bin containing the rigidity of the event. An MDR efficiency is 
thus calculated for each sub-bin with a center rigidity Rj and then a mean efficiency 
in bin i (gj) is determined by weighting the efficiency in the sub-bin j (sj) by a 
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theoretical flux, which is ~ R~ 2 7 , in the following way: 

^ £?=i RJ 2 - 7 x 6J x A Rj ^ ^ 



As shown in figure [47T8] the MDR efficiency decreases dramatically above ~20 GV. 




Figure 4.18. The MDR efficiency for protons evaiuated from flight data. The blue 
points show the efficiencies in sub bins. The red points show the mean efficiency in 
wider bins. 



4.2.7 Trigger efficiency 

The PAMELA trigger efficiency is a product of the trigger efficiency for each ToF 
layer required in a particular trigger configuration. The trigger efficiency is calcu- 
lated to exceed 0.997 with an error of the order 0.5 x 10 -4 [140] . Compared to other 
factors discussed in this section, the trigger inefficiency is completely negligible and 
therefore will be omitted when calculating the total selection efficiency. 



4.2.8 Total selection efficiency 



All grouped efficiencies discussed above are compared in figure 4.19 The tracker 
efficiency here is defined as the total efficiency of the basic and additional tracker 
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selection cuts, calculated by multiplying the basic and the additional tracker effi- 
ciencies. The selection efficiency is dominated by the tracker selection. Above 100 
GV, the MDR selection also play an crucial role. After all individual efficiencies are 
derived and possible correlations are understood, the total antiproton and proton 
selection efficiencies can thus be calculated by multiplying all terms and associated 



errors, as shown in figure 4.19 
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Figure 4.19. Top: All grouped efficiencies integrated over the whole live time. Bot- 
tom: Total efficiency integrated over the whole live time for protons and antiprotons. 
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4.3 Correction factors 

Apart from selection efficiencies, several other correction factors, i.e. the geometri- 
cal factor and hadronic interaction losses in the instrument, the live time of mea- 
surements and transmission through the geomagnetic field, should be applied to 
the selected candidates to reconstruct the antiproton flux. Moreover, to derive 
the antiproton-to-proton flux ratio, differences between the corrections between 
antiprotons and protons are also addressed. 



4.3.1 Geometrical factor and hadronic interaction losses 
Geometrical constraint 

The flux intensity is generally a product of a count rate and a proportionality factor 
called the gathering power of the detector. Under the hypothesis of an isotropic 
flux, the gathering power is expressed as a geometrical factor , defined by 



G(R)= / dCl / dS\cos9\f(x,y,9,4>,R), (4.8) 
Jn Js 

where R is the rigidity, f2 is the total solid angle, S is a reference plane orthogonal 



to the z axis defined in figure 3.2 and / is a weighting function that is either f 
or depending whether the trajectory of incident particle satisfies the acceptance 
requirements of the apparatus. The acceptance requirements are: 

• the trajectory must cross at least one of the two layers in each plane of the 
ToF system: (Sff OR S12) AND (S21 OR S22) AND (S31 OR S32). 

• the particle must cross all the 6 planes of the tracking system. 

• the trajectory must be fully contained in the magnetic cavity without touching 
the walls of the cavity. 

The geometrical factor is dependent on the rigidity of the incident particle. As 
lower rigidity particles are more deflected by the magnetic field towards the walls 
of the magnetic cavity, where they are absorbed before reaching the lower face of 
magnetic cavity, the geometrical factor is expected to decrease at low rigidities. 
At high rigidities where particle trajectories are approximately straight, the geo- 
metrical factor is expected to be constant. The geometrical factor presents the 
geometrical constraints of a particle telescope and docs not depend on the particle 
species. 

In order to calculate the geometrical factor, an approach based on the work by 
Sullivan [146] has been performed with simulations. A set of particles are generated 
on a generation surface just above the scintillator plane SI, each with random 
parameters {x, y, 9, (ft, i?), in which (x, y) is the coordinate, the (6, </>) is the incident 
direction and R is the initial rigidity. 
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The measured PAMELA magnetic field in the spectrometer has been used in 
the simulation. For a given rigidity, the geometrical factor is given by: 

Gf = — — ■ Gg erl , (4.9) 
ntot 

where n se i is the number of selected particles satisfying all the acceptance require- 
ments, ntot is the total number of generated particles, and G gen is the gathering 
power of the generation surface with area S gen . G gen can be expressed as: 

Gg en = I dCl dS\cos0\ = S gen n (l - cos 2 9 max ) , (4.10) 

where 6 ma x is the maximum generation zenith angle. As only downward-going 
particles are of interest, the angular domain is limited to the downward hemisphere, 
characterized by7r/2<#<7r. 

Interaction losses 

The geometrical factor discussed above depends only on the geometrical constraints 
of the instrument. However, due to the presence of material a particle traversing the 
apparatus may interact inside the acceptance. Since an incident particle entering 
the acceptance of PAMELA should always be accounted for in the calculation of 
the particle flux, inelastically interacting events, which are rejected by the selection 
criteria applied on flight data should be evaluated. In detail, two effects should be 
considered: 

• the loss of particles which would traverse the acceptance cleanly if not inter- 
acting inside the acceptance. 

• the gain of particles which would not traverse the acceptance if not scattering 
inside the acceptance. 

This correction is not accounted for in the selection efficiency calculation since 
particles which interact inside the instrument and produce secondaries are rejected. 
Therefore, the geometry of the apparatus and all physical processes are implemented 
in the simulation to estimate the effective acceptance which includes the correction 
due to inelastic interaction effect. An independent work has been performed by 
Bruno [149] , A sample of downgoing particles have been isotropically generated 
from a surface placed above the "dome" (the top container of the instrument) . The 
starting point, direction and the area of generation surface were chosen to ensure 
no bias in the calculation of the fraction of in-acceptance events and include those 
events in the sample whose tracks were initially not contained in the acceptance 
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but were deflected into the acceptance by scattering (mainly in the dome). The 
effective acceptance is calculated by: 

A F = ^L. Ggen . (4.11) 
n-tot 

Since the interaction processes are included in the simulation, the particle loss and 
gain due to interactions are naturally included in the effective acceptance. 
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Figure 4.20. The PAMELA effective acceptance for protons (blue) and antiprotons 
(red). From 149 . 



The results of the calculation are shown in figure 4.20 both for positively- 



charged particles and negatively-charged particles |149j . Below 1 GeV, the geome- 
try significantly constrains the acceptance. However, above that energy, the shape 
of the effective acceptance mainly reflects the cross sections of inelastic interactions 
for protons and antiprotons. The inelastic cross section for protons increases with 
energy from 1 GeV to 2 GeV and remains constant between 2 GeV and 1000 GeV, 
resulting a decreasing effective acceptance and an almost unchanged acceptance 
below and above 2 GeV. For antiprotons, the inelastic cross sections increases dra- 
matically until 100 GeV, resulting in an increasing effective acceptance below that 
energy. The difference between protons and antiprotons is due to the difference in 
inelastic cross section for these two species. The hadronic generator FLUKA is em- 
ployed in GPAMELA simulation code |147j . which was developed by the PAMELA 
collaboration based on the GEANT package [148] version 3.21. In order to simulate 
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the hadronic physical processes and to generate secondary cascades from hadron- 
nucleus interactions, different models are used in FLUKA generator depending on 
the energy of the cosmic ray hadrons. An error of 5% is added to the effective 
acceptance considering systematic uncertainties related to the hadronic generator. 

4.3.2 Transmission through the geomagnetic field and live 
time 



As discussed in section 4.1.5[ only particles satisfying the following requirement are 
accepted: 

rigbin lowerlimit > cutoff pamela = 1.3 x SVC. 

Galactic particles not satisfying the rigidity requirement arc rejected and should be 
accounted for when deriving the flux. 

This correction was calculated together with the live time, which is defined as 
the time when the experiment is operational and ready for a new trigger. Contrary 
to the live time, the time when the instrument is switched off or is reading out and 
processing data is called the dead time. The live time corrected by transmission 
through the geomagnetic field is calculated by: 

rbini ower ii m u/1.3 

T Uve (bin)= / f(SVC)dR, (4.12) 

Jo 

where Tu ve (bin) represents the live time spent for 1.3- SVC lower than bini ower u m i t 
(the low edge of the bin), and / (SVC) is the cutoff distribution weighted by relative 
live time. Assuming that the particle flux is isotropic, the loss of particles can be 
compensated for by multiplying the measured particle intensity with the inverse of 



the corrected live time. The result is presented in figure 4.21 Below about 20 GV 



a continuous increase of live time with rigidity can be seen. 
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Figure 4.21. The live time spent for PAMELA cutoff < bini ower i imit . 

4.4 Antiproton flux and antiproton-to-proton 
flux ratio 

The raw number of antiprotons and proton candidates surviving all the selection 
criteria are presented in in table |4.1| These candidates are selected by accounting 
for the time when the tracker was not operating satisfactorily or the statistics is 
too low, as discussed in section |4.2.1| The days of detector failure are labelled 
as "non-operational" days and are excluded in the data analysis to ensure stable 
efficiencies. For the rigidity range considered in this thesis, i.e. between 2.23 GV/c 
and 180 GV/c, the pion contamination was estimated to be negligible |149j . as 
well as the electron contamination |141j . Possible proton contamination due to the 
spillover effect was studied using both simulation and flight data. The antiproton 
selection criteria except the MDR cut is applied to simulated protons to reproduce 
the spillover observed in real flight data, and the proton contamination after ap- 
plying the MDR cut is then determined. The systematic uncertainty due to the 
proton contamination only exists in the highest energy bins and was estimated to 
be -20% for the rigidity bin 48.5-100 GV/c and -30% for the bin 100-180 GV/c. 

In order to construct the antiproton flux at the top of the PAMELA payload, 
the raw number should be corrected for selection efficiencies, hadronic interaction 
losses and the geometrical factor, transmission through the geomagnetic field and 
measurement live time. The selection efficiencies discussed in section 14.21 show a 
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Rigidity GV/c 


Antiprotons 


protons 


2.23 - 


2.58 


49 


1701576 


2.58 - 


2.99 


77 


1644062 


2.99 - 


3.45 


78 


1485846 


3.45 - 


3.99 


96 


1376729 


3.99 - 


4.62 


102 


1246995 


4.62 - 


5.36 


107 


1130713 


5.36 - 


6.23 


108 


993868 


6.23 - 


7.27 


98 


873283 


7.27 - 


8.53 


85 


762115 


8.53 - 


10.1 


94 


671588 


10.1 - 


12.0 


105 


567089 


12.0 - 


14.6 


78 


527378 


14.6 - 


18.1 


64 


359288 


18.1 - 


23.3 


55 


337924 


23.3 - 


31.7 


41 


272516 


31.7 - 


48.5 


36 


194277 


48.5 - 


100.0 


22 


106888 


100.0 - 


180.0 


3 


13596 



Table 4.1. The antiproton and proton candidates after excluding "non-operational 
days. 



time dependence, especially the tracker efficiency. For a bin i, if the raw number 
in each day Nj is corrected with a daily efficiency, ej , the systematic effects due to 
the time variation of detector performance can be minimized. However, as there 
are a small number of antiproton candidates per day, if an event is selected by 
chance in a day when the efficiency is low and thus no candidate is expected to be 
selected, the corrected number will be overestimated. Therefore, instead of using a 
daily efficiency, the correction is done on the basis of a 'bunch' of several days to 
eliminate this overestimation by increasing the statistics, as well as to reduce the 
effect of time dependent efficiencies. A value of 60 days is chosen to be appropriate 
for one bunch in final calculation. 



The resulting flux is presented in table 4.2 The range of an energy bin is 
converted from the range of rigidity bin in table |4.1| assuming a singly-charged 
particle and by using Ek = \fp 2 c 2 + m 2 c 4 — mc 2 , where Ek is the kinetic energy, m 
the mass of the particle, p the momentum and c the velocity of light. Data points 
are centered in each bin according to a technique developed by Lafferty and Wyatt 
(150] . For a bin with lower limit Ei and upper limit E u , the weighted center E c is 
determined as the abscissa value at which the measured spectrum is equal to the 
expectation average value of the "true" spectrum, which can be expressed as 

1 ^ 



f( E J=w — w f( E ) dE > ( 4 - 13 ) 

- &l J Ei 



<S(i 
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where / (E) is a theoretical model of antiproton spectrum |151| 1103] . 

Since for antiprotons and protons the discrepancy between correction factors 
only exists for the effective acceptance and the calorimeter efficiency due to the 
difference of hadronic cross sections, the p/p ratio can be calculated as 



R {bin) 



N p (bin) I (epcaio (bin) x G p ) 
N p (bin) I (Ep Ca i (bin) x G p ) ' 



(4.14) 



where Np is the total number of antiproton candidates, N p the total number of pro- 
ton candidates, £p ca io the calorimeter efficiency for antiprotons, £ pcQ z the calorime- 
ter efficiency for protons, Gp the acceptance for antiprotons and G p the acceptance 
for protons. The resulting ratio, which increases with energy up to ~ 10 GeV and 
then flattens, is presented in table [472] 
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Table 4.2. The antiproton flux and p/p flux ratio measured by PAMELA. The 
first set of errors refer to the la statistical errors and the second set are systematic 
errors. 



The statistical and systematic errors are shown separately in table |4.2| For both 
the antiproton flux and the p/p ratio, statistical errors dominate over the entire 
energy range. Figure 4.22 and figure 4.23 compare the antiproton flux and the p/p 
ratio derived in this work respectively with other experimental data. PAMELA 
measurements are consistent with other measurements but with significantly better 
statistics. The results derived in this work are also compared with those officially 
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published by the PAMELA Collaboration |152j extending to lower energies. The 
agreement is excellent^] 

The majority of cosmic ray antiprotons are generally believed to be produced 
secondarily by interactions of cosmic ray protons with the ISM. Due to the kine- 
matic constraints on the antiproton production, a peak around 2 GeV is expected to 
appear in the antiproton flux and is confirmed by results shown here. The antipro- 
ton spectrum is a useful tool to constrain propagation models. Different models 
may give different predictions on antiproton fluxes. For example, the diffusion 
reacceleration models usually produce too few antiprotons [HI [75] . Some models 
expect slightly different antiproton spectra which could not be discriminated by 
antiproton data published before PAMELA (e.g. as in |103j ). Thanks to the more 
accurate antiproton data provided by PAMELA, stronger constraints may be able 
to be placed on cosmic ray propagation models. In the next chapter, the antiproton 
flux and the p/p ratio measured by PAMELA will be used to further study cosmic 
ray propagation models. 



2 The published results used a different method to estimate the total efficiency. Loose selections 
were also used in the rigidity range 6.23-14.6 GV once it became clear that pion contamination 
was minimal. 
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Figure 4.22. The p flux measured by PAMELA compared with the results from 
other experiments (see 86 (and references therein) 11531 and |154l ). PAMELA mea- 
surements related to this work are shown in red symbols. PAMELA measurements 
published in |152| are shown in blue symbols. 
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Figure 4.23. The p/p ratio measured by PAMELA compared with the results 
from other experiments (see references in 141 ). PAMELA measurements related to 
this work are shown in red symbols. PAMELA measurements published in I152| are 
shown in blue symbols. 
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4.5 Discussion 

Before the PAMELA experiment was carried out, the CAPRICE1998 measurement 
indicated an increasing trend in the antiproton flux above 10 GeV and therefore 
a possible presence of a primary antiproton contribution. The much more precise 
PAMELA results disfavor this trend and generally agree with the pure secondary 
production calculations. However, in contrast to antiproton flux, PAMELA ob- 
served a clear rise above 10 GeV in the positron fraction disagreeing with the 
prediction of secondary positron production. A considerable number of interpre- 
tations, including astrophysical or exotic primary sources, have been proposed to 
explain both the antiproton data and the unexpected positron fraction. Nearby 
Pulsars have been proposed as a good candidate for the positron excess since 
electron-positrons pairs are expected to be produced in pulsars but no antipro- 
tons (e.g. |155[I156(I157| ). Another explanation suggested for the steep increase in 
the positron fraction is provided by positrons and electrons as secondary products 
of hadronic interactions inside aged SNRs |158j . This model also predicts an an- 
tiproton flux compatible with the PAMELA data but predicts a harder antiproton 
component beyond the 100 GeV energy region which strongly differs from the con- 
ventional antiproton production created only from spallation in the ISM [159 . The 
dark matter scenario, as discussed in section |1.3.2| may also result in anomalous 
features in spectrum of cosmic ray antiparticles. Both antiprotons and positrons 
are generally believed to be the final states of dark matter annihilation, however, 
no excess over the standard secondary background is seen in the antiproton spec- 
trum and this has therefore placed strong constraints on dark matter properties. 
The dominant final states of dark matter annihilation or decay are leptons instead 
of quarks [551 U60| . Candidates like the Kaluza-Klein particles which annihilate 
mostly to leptons, provide the possibility to reproduce the observed electron and 
positron spectrum [S3 [M] . Dark matter annihilation channels of charged gauge 
bosons (W, Z) can only be accommodated in PAMELA data with > 10 TeV dark 
matter mass [160] , unless electrons and antiprotons have significant different boost 
factors (for example in [161] ). Nevertheless, a robust estimation of the cosmic ray 
propagation is a foundation of the investigation on possible primary sources and 
will be studied in next chapter. 



Chapter 5 

Constraints on transport and 
acceleration models 



As discussed in the last chapter, the cosmic ray antiprotons measured by PAMELA 
can be considered as an important means to test propagation models and dark 
matter properties. In previous works published in the literature (see references 
in section 2.5), cosmic ray propagation is usually studied using data from differ- 



ent experiments. However, inconsistencies might exist between data sets and thus 
introduce systematic errors in the final results. Since PAMELA provides measure- 
ments of a variety of cosmic ray species, this systematic effect can be potentially 
reduced. The work presented in this chapter examines whether PAMELA antipro- 
ton and proton data can provide strong and reliable constraints on the source and 
propagation parameters. In addition, the value of the upcoming PAMELA B/C 
ratio is demonstrated by employing the B/C ratio from other experiments cover- 
ing a comparable energy range to PAMELA. The unprecedented accuracy of the 
PAMELA B/C ratio is expected to give even better constraints. 

The propagation and source parameters under study are summarized in sec- 
tion 15.11 The data used in this work are shown in sectio n l5~2l Two statistical 
approaches, i.e. the \ 2 minimization method (see section 5.3) and the Bayesian 



method (see section 5.4), are used to constrain the propagation and acceleration 
models. Moreover, the Bayesian analysis allows us to test hypotheses, for which the 
electron flux and the positron fraction are calculated as a consistency check and as 
an input for future dark matter searches. 



5.1 Summary of studied parameters 

Different propagation models are studied in this work (as also defined in chapter 
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• the plain diffusion models (PD); 

• the diffusion reacceleration models (DR); 

• the diffusion convection models (DC); 

• the diffusion reacceleration convection models (DRC). 

The public numerical package GALPROP (version 54.0.572) |117] is used. In 
GALPROP, the source abundance of protons is normalized based on the propagated 
proton spectrum at solar position for an energy of 100 GeV, referred to as N p . The 
normalizations of other nuclei are then scaled by their source abundances relative 
to that for protons. Therefore, instead of the absolute normalization abundances of 
the injection spectra for different cosmic ray species, N p is used as a free parameter 
to characterize the source term. Other free parameters related to the source term 
and the propagation processes were described in section |2.3| A summary of these 
model parameters follows here: 

• Dq: a free normalization of the diffusion coefficient at a reference rigidity 
4 GV; 

• (5: the spectral index of the diffusion coefficient; 

• va- the Alfven speed characterizing the reacceleration effect; 

• dV/dz: the derivative of convection velocity; 

• v. the injection spectrum index; 

• N p : the normalization of the propagated proton spectrum at a reference ki- 
netic energy of 100 GeV. 

For stable cosmic ray species, a degeneracy exists between Do and the halo size 
of the Galaxy, Zh- Therefore, these two parameters cannot be constrained simul- 
taneously by using stable secondary-to-primary ratios. In this work, Zh = 4 kpc is 
assumed in agreement with earlier GALPROP studies of 10 Be, 26 Al, 36 C1, 54 Mn, 
and the B/C ratio [7Tlll70j and to ease the comparison of results. Other GALPROP 
parameters are held at the conventional reacceleration configuration [7Tlll03j . tuned 
to reproduce the ACE isotopic abundances of |171j . For studies of the B/C ratio, 
the nuclear chain starts from 28 Si since all primary elements from Si down to C 
have an important effect on the B/C ratio. For studies only including proton and 
antiproton data, the nuclear chain starts from He since primary helium nuclei 
and proton interactions with the ISM are dominant in the secondary production of 
protons and antiprotons. 

To reproduce solar modulation, the force-field approximation described in sec- 
tion |l.l| is used in this analysis. The modulation potential $ is chosen to follow that 
reported by each experiment, i.e. 500 MV for PAMELA [162] and HEA03 [59], 
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325 MV for ACE-CRIS [165], 850 MV for CREAM- 1 [166], 400 MV for Spacelab-2 
[55] and 450 MV for AMSOlfTBT]. It is worth to note that the derived values of $ 
suffer from some uncertainties since their determination is phenomenological and 
depends on the choices of interstellar spectra (see equation |1.1[ ) . To take into ac- 
count these uncertainties, the solar parameters are free in the Bayesian analysis. 
The prior probability distributions (simply called priors), stating our initial knowl- 
edge on the parameters before seeing the data, are specified as described in section 
|5.4.1| Unlike the model parameters, the solar parameters can be considered as nui- 
sance parameters which are include in the parameter scan but are not of primary 
interest. 



5.2 Data 

As mentioned, PAMELA achieves significantly better statistics and extends the 
energy range compared to previous experiments, especially on antiparticles, e.g., 



the antiproton flux and the p/p ratio, as shown in figure 4.22 and 4.23 These 
measurements cover the energy range from 60 MeV to 180 GeV. Such precise mea- 
surements enable us to put better constraints on cosmic ray transport parameters. 
In the analyses presented in this chapter, the published PAMELA data on the p/p 
ratio [152 are used to be consistent with other studies in the literature. Together 
with the antiproton data, we also use the proton flux measured by PAMELA with 



great accuracy from 400 MeV to 1.2 TeV |162) . as shown in figure 5.1 to improve 
constraints on the primary injection spectrum. A hardening in the spectra around 
200 GeV can be seen. Some interpretations were proposed to explain this hard- 
ening, for example dispersion in the source injection spectra [168 or the neutral 
atoms presented during the acceleration process in a shock |169j . This spectrum 
hardening, however, is not modeled in this work in order to focus on the propagation 
processes described in chapter [2] 

We propose to exclusively use PAMELA data to study cosmic ray propagation 
for several reasons. Usually in order to constrain propagation parameters, it is nec- 
essary to combine data sets from a variety of experiments to cover a wide enough 
energy range. As pointed out in [73] . errors might be underestimated for an exper- 
iment. In order to compensate the systematic discrepancies in the reported uncer- 
tainties from data sets, one can therefore introduce a set of nuisance parameters 
to rescale the reported errors. These rescaling factors increase the computational 
time for parameter space scans and hinder model selection. These difficulties can 
be avoided by using only PAMELA data, which additionally allows the use of a 
smaller parameter set. Another unavoidable problem caused by incorporating data 
sets from various experiments is that the modulation potential <3 ) based on the as- 
sumed interstellar spectrum of cosmic ray species may differ between experiments. 
Relatively poorly understood solar physics makes studies of cosmic ray transport 
in the Galaxy more difficult. Using only PAMELA data decreases uncertainties on 
derived propagation parameters by including 'I'pamela as a nuisance parameter, 
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Kinetic energy (MeV) 



Figure 5.1. Proton flux measured by PAMELA above 1 GeV/n 11621 . 



though a potential bias might still arise from the simplified approximation of solar 
modulation used. 

The B/C ratio, one of the quantities most sensitive to the propagation parame- 
ters, is expected to be measured from 100 McV/n to 200 GeV/n by the PAMELA 
experiment in near future. PAMELA is also able to measure hydrogen and helium 
isotopes, i.e. 2 H and 3 He, over an energy range from 100 MeV/n to 700 MeV/n 
and 900 MeV/n respectively [1631 1164] . These isotopes are believed to be pro- 
duced during interactions of 4 He. The ratios 2 H/ 4 He and 3 He/ 4 He have been 
shown to be as constraining as the B/C ratio |175j . However, all these ratios from 
PAMELA are not available yet, therefore the B/C ratio measured by previous ex- 



periments [551 1591 HB31 11661 1167| (see figure 5.2 1 are used in this work to constrain 
the transport parameters. The energy range covered by the B/C data sets chosen 
in the analysis is comparable to the energy range that PAMELA will provide. The 
ratios of 2 H/ 4 He and 3 He/ 4 He are not included in the analyses for three reasons. 
Firstly, they are less accurately measured by individual experiment than the B/C 
ratio because of the detetors' isotopic separation ability. Secondly, including more 
data from various experiments will increase the uncertainties due to discrepancies 
between data sets. Thirdly, a larger number of modulation parameters need to be 
dealt with. 

Different combinations of data sets are used in the analysis presented in this 
chapter and are labeled as follows: 
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Figure 5.2. B/C ratios measured by HEA03 [59], ACE-CRIS [TBS], CREAM-1 
[1551 . Spacelab-2 [58] and AMS01 [T571 . PAMELA is expected to measure the B/C 
ratio from 100 MeV/n to 200 GeV/n (the red band). 



"-o" : only the antiproton flux; 

"-a": only the B/C ratio; 

"-b": a combination of the p/p ratio and the proton spectrum; 

"-c": a combination of the B/C ratio, the p/p ratio and the proton spectrum. 

5.3 x 2 minimization approach 

To understand how well a model reflects the observed data, a commonly performed 
statistical test is the % 2 'goodness of fit' test. Assuming we have a number of N 
data points (JQ, Yi) with Gaussian distributed errors, the x 2 test statistic is defined 
as: 



x 2 (©) 



^ (/(x a ,0)-y t ) 2 
»=i 



07 



(5.1) 



where / (Xj, ©) is the theoretical value at abscissa Xj, and <7j are the uncertainties 
on the measurements Y^. Since the theoretical expectation depends on the vector 
of physical parameters ©, the best-fit parameters can be extracted by minimizing 
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the x 2 . Usually the reduced x 2 is used, i.e. x 2 divided by the number of degrees 
of freedom (d.o.f), where the d.o.f is equal to the number of data points minus the 
number of free parameters. 

GALPROP is interfaced with the minimization library MINUIT [172] for this 
analysis. The MINUIT package provides several algorithms to perform a mini- 
mization of a multi-parameter function, such as the the Nelder-Mead SIMPLEX 
algorithm |173j and the MIGRAD algorithm based on the variable metric method 
|174j . An error matrix is calculated as a by-product of MIGRAD and the parabolic 
errors are derived. More general asymmetric errors can be further produced by 
the MINOS method. In this work, / (A;, 0) corresponds to the cosmic ray fluxes 
and/or flux ratios given by GALPROP at kinetic energy Xi. The \ 2 function is 
derived from the GALPROP prediction of / pQ, 0) and the experimental data 
described in section 15.21 The efficient MIGRAD method is used to evaluate the 
best-fit parameter values and MINOS is used to derive the uncertainties of the 
parameters. If MIGRAD fails to converge then the minimizer is switched to the 
slower SIMPLEX method. 

5.3.1 Analysis and results for an unmodified diffusion 
coefficient 

The antiproton data measured by PAMELA and the B/C ratio reported by other 
experiments are used separately to check their constraining power for propagation 
parameters. Instead of the p/p ratio, the antiproton flux is used since once the 
proton flux is fixed, the antiproton flux is more sensitive than the p/p ratio to the 
propagation parameters. The source injection index is fixed at v — 2.3 and the prop- 
agated proton spectrum is normalized as N p = 4.69 x 10 -9 cm -2 sr _1 s _1 MeV -1 
at 100 GeV to fit the PAMELA proton data. The simplest model, i.e. PD, 
is studied. Results obtained from the antiproton flux (PD-o model) are D = 
3.88±0.14 x 10 28 cm 2 /s and 5 = 0.479± 0.024, which are close but less constraining 
than those obtained from the B/C ratio (PD-a model), D = 5.57±0.04x 10 28 cm 2 /s 
and 5 — 0.490 ±0.008. By fitting the antiproton flux, the best-fit PD-o model gives 
a reduced x 2 of 0.70, which means that the data are overfitted by the PD-o model. 
The DR-o model, with = 14±^2 kms -1 and diffusion coefficients consistent 
with the PD-o model, does not change the value of \ 2 - The large error on indi- 
cates that reacceleration is not constrained by only using the measured antiproton 
spectrum. The convection velocity dV/dz is converged at zero, since the model 
with convection (DC-o) always increases the x 2 value compared to the PD-o model 
and is disfavored. The PD-o model is sufficient here to reproduce the antiproton 
data due to two possible reasons: (1) the low energy antiprotons are primarily pro- 
duced through "tertiary" processes coming from inelastic scattering of high energy 
antiprotons, therefore the antiproton flux may not be very sensitive to the other 
low energy processes, i.e. reacceleration and convection; (2) the antiproton flux is 
dependent on the injection spectrum of primaries since antiproton production is 
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significantly influenced by the kinetic energy of their progenitors, i.e. protons and 
helium nuclei. 

Tighter constraints are placed on transport parameters by using the B/C ratio, 



which are summarized in table 5.1 as well as the corresponding x 2 /d.o.f. From 
the reduced % 2 i PD-a and DC-a models can not reproduce the B/C ratio. Large 
deviations at low energy can be seen in figure |5.3| from the comparison of the data 
with the best-fit PD-a and DC-a models. The reacceleration process (DR and DRC 
models) can explain the B/C ratio adequately. The DRC-a model is found to best 
fit the B/C data. Compatible results are obtained by fixing the injection spectrum 
index v = 2.5, which indicates that the sensitivity of the B/C ratio to the injection 
spectrum is weak. This is expected since the boron nuclei have almost the same 
energy /nucleon as their primaries. 




Figure 5.3. The B/C ratio for the best-fit parameters of PD-a, DR-a, DC-a and 
DRC-a models as listed in table I5TT1 

In order to obtain complementary information on transport and source param- 
eters, a simultaneous fit to both the secondary-to-primary ratios and the primary 
fluxes is necessary. Using only PAMELA data, the proton flux is combined with the 
p/p ratio to estimate the parameters. When only the data of antiproton spectrum 
were fitted, the source parameters were fixed and the simplest model, PD-o, seems 
sufficient to describe the antiproton data. However, if the source parameters are 
varied, the derived values of transport parameters will be considerably changed. 
Therefore a simultaneous fit on the p/p ratio and the proton flux may allow us to 
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Figure 5.4. The p/p ratio (top) and the proton spectrum (bottom) for the best-fit 
parameters of PD-b and DC-b models as listed in table |5.f | 
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constrain all the processes. The results listed in table |5.1] show that parameters for 
both PD-b and DC-b models can be constrained. However, the reliability of these 
best-fit parameters needs to be questioned. A rather high spectral index of the dif- 
fusion coefficient S = 0.84 ± 0.04 appears in the PD-b model, which can reproduce 



the proton spectrum but not the p/p ratio, as shown in figure 5.4 This bias on 



the estimated parameters is due to the dominant weight of the proton flux since it 
is more precisely measured than the p/p ratio. The DC-b model, with the reduced 
X 2 close to 1, gives comparable best-fit parameters to that for the DC-a model and 
can generally fit both the proton flux and the p/p ratio. Contrary to the indication 
from the B/C ratio, reacceleration is disfavored by fitting simultaneously the p/p 
ratio and the proton flux, i.e. va — > for DR-b and DRC-b models. Generally, the 
derived best-fit propagation parameters have errors at least twice larger than the 
ones obtained from the B/C ratio, showing that the combination of p/p ratio and 
the proton flux is not as constraining as the B/C ratio. 

Since the estimated parameters might not be reliable by fitting only the p/p ratio 
and the proton flux, the B/C ratio is also included in the simultaneous fit. The B/C 
ratio is more constraining than the p/p ratio and therefore will decrease the bias 
on the transport parameters. The best-fit PD-c and DC-c models can reproduce 
all the data except for the low energy B/C ratio, as shown in figure 5.5 Unlike the 
PD-b model, the PD-c model gives a lower spectral index of the diffusion coefficient 
S = 0.495 ±0.007 which provides a satisfactory fit on the p/p ratio. Reacceleration, 
which is expected to explain the low energy B/C ratio, is still disfavored since it 
conflicts with the proton flux. One way to solve this disagreement is to introduce 
an unphysical adhoc break in the injection spectrum, as did in |44l 173] . The same 
configuration, referred to as "DR II", is also tested here to fit the B/C ratio, the 
p/p ratio and the proton spectrum. The best-fit parameters of this DR II-c model 



are also given in table 5.1 However, the predicted B/C ratio of this model is still 



higher than the data below I GeV (see figure 5.5). This is because that a higher 
Va (as the best-fit value of the DR-a model) than the one obtained for the DR II-c 
model, is needed to explain the B/C data. In order to account for the low energy 



B/C ratio, a nonlinear diffusion coefficient will be further studied in section 5.3.2 



Comparison with previous studies 

Except for the biased values obtained for PD-b and DR-b, the spectral index 5 
of the diffusion coefficient is well constrained between 0.3 and 0.65 for all models 
considered. Whereas the PD models favour a Kraichnan turbulence spectrum of 8 = 
0.5, the DC models favour a slightly higher value of 5 between 0.62 and 0.65. The 
Kolmogorov spectrum of turbulence, 5 = 1/3, is only recovered for the DR-a model 
and DR II-c, in agreement with earlier studies (e.g., [3H[73]). Including observations 
of primary nuclei, i.e. the proton flux, tends to disfavour reacceleration unless a 
break is introduced in the injection spectrum as done in the DR II-c model. The 
Galactic wind velocities obtained for the DC models are around 10 kms _1 kpc _1 
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and in good agreement with other studies, e.g., 44 . The same is valid for the 
Alfven velocities of about 35-40 kins -1 for the DR-a and DRC-a models. 



5.3.2 Analysis and results for a modified diffusion 
coefficient 

The discrepancy between the low energy B/C observations and the models which 
are compatible with other data, could be explained by nonlinear MHD wave effects, 
i.e. the spectral wave density decreases from small to large wave numbers [1031 . 
This turbulence dissipation effect can result in a low energy dependence of the 



cosmic ray diffusion coefficient, as mentioned in section 2.5 A spatial diffusion 
coefficient D xx = Dof3 v {ft/ po) is adopted to account for this effect, where r\ is 
added as a free parameter in the fit. Unlike introducing an artifical break on the 
diffusion coefficient at a rigidity 3 GV or 4 GV as used in [531131], or on the injection 
spectrum at 10 GV used in [HI [73], this approach is physically motivated. 

Models including rj to parameterize low energy MHD physics are referred as PD, 
DR, DC and DRC. The same combinations of data sets as used in the previous 
section, are employed to study these models. All the results are summarized in 
table |5.2| The effect of r\ competes with effects of convection and reacceleration. 
If r) is included in the fit of the B/C ratio, no convection is favored. To fit both 
the p/p ratio and the proton flux, the values of the reduced x 2 are much smaller 
than 1.0. This unexpectedly good fit may be caused by the possible correlated 
and/or overestimated errors, or the errors are not Gaussian distributed, which then 
would result in a significant deviation from the % 2 distribution and make the \ 2 
minimization approach inappropriate in this specific case. A more precise study is 
beyond the scope of this thesis, but will be addressed in a future work. The smaller 
values of the reduced % 2 indicates that models including a low energy dependence 
in the diffusion coefficient are preferred to explain the p/p ratio and the proton 
flux. But using only the p/p ratio and the proton flux is not enough to constrain 
models including rj. In the simultaneous fit of the B/C ratio, the p/p ratio and the 
proton flux, the convection velocity is converged to zero due to degeneracy between 
rj and dV/dz and the dominant effect of rj at low energy. 

The PD-c model is generally consistent with all the data except for a slight 
overprediction of the p/p ratio below 10 GeV. The best-fit values of S = 0.621±0.009 
and rj = — 1.75±0.10 are close to the values given in |103j . i.e. S = 0.60 and rj = —2. 
Reacceleration can be invoked to fit the data. However, since the required Va is 
very weak and the reduced x 2 for model DR-c is compatible with the one for model 
PD-c, reacceleration seems not to be necessary to explain the data. This can also be 



illustrated in figure 5.6 Nevertheless, 77 is found to dominate over other competing 
processes and may result in too many degenerated parameters being constrained 
at low energy. 
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5.3.3 What did we learn from the % 2 study? 

From the % 2 study, the following conclusions can be drawn: 

• The antiproton data alone is not enough to constrain different propagation 
processes. 

• Reacceleration can explain the B/C data but produces too many protons at 
energies of a few GeV. A break in the injection spectrum is included in some 
studies to fit the data of the proton spectrum. If this break is not introduced, 
however, in order to fit all the data, including the B/C ratio, the p/p ratio 
and the proton spectrum, reacceleration is always disfavored, i.e. va — > 0. 

• Except for the low energy B/C ratio (< 1 GcV/n) reported by ACE-CRIS, 
the PD and DC models can generally explain the high energy B/C ratio, the 
p/p ratio and the proton spectrum. 

• A low energy dependence applied to the diffusion coefficient can fit all the 
data. But other important processes at low energy, such as convection and 
reacceleration, are not allowed to be studied since the low energy dependence 
of the diffusion coefficient dominates over these effects at low energy. 

• The estimated parameters might be biased in the simultaneous fit to both the 
secondary-to-primary ratios and the primary fluxes due to the very precise 
proton data. The parameters characterizing the low energy precesses could 
also be biased due to the solar modulation which is simply modeled by an 
effective parameter, used in the force-field approximation. 

• The statistical uncertainties on parameters depend on the parameters under 
study and the data used. The errors on N p are constant to be about 1% since 
this parameter is normalized to the measured proton flux. By fitting all the 
data simultaneously, parameters other than N p have errors less than 10%, 
which are at least twice precise than those estimated by fitting PAMELA 
proton and p/p data. However, if the parameter rj is included in the fit, the 
values of va and dV/dz are not constrained very well, i.e. uncertainties are 
generally larger than 25%. 
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Figure 5.5. The B/C ratio (top), the p/p ratio (middle) and the proton spectrum 
(bottom) for the best-fit parameters of PD-c, DC-c and DR II-c models as listed in 
table HI] 
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Figure 5.6. The B/C ratio (top), the p/p ratio (middle) and the proton spect rum 
(bottom) for the best-fit parameters of PD-c and DR-c models as listed in table 5.2 
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5.4 Bayesian approach 

As discussed in the last section, a simultaneous fit to both primary and secondary 
cosmic ray data might bias the parameters. This was also argued in |175j . Indeed, 
primary fluxes are more prone to systematics and are more sensitive to solar mod- 
ulation than secondary-to-primary ratios. This bias can be reduced by specifying 
priors on the source parameters and taking into account the uncertainties on the 
solar modulation potentials. For these reasons a Bayesian method is used, allowing 
parameters to be estimated based on prior knowledge and information contained 
in the likelihood, i.e. the probability to observe the data measured for a particular 
model assumption. 

Given the observed data set D and the parameters © under study in a hypothesis 
(model) H, Bayes' theorem states that: 

P(D\H) ' V ' 

where P(@\T),H) is the posterior probability density function (p.d.f.) of the pa- 
rameters, P (D|0, H)=L (0) is the likelihood, P (&\H) is the prior, and P (D\H) 
is the Bayesian evidence. 

The posterior probability distributions of the transport and source parameters 
are derived by Bayesian inference which can naturally produce credible regions 
in the parameter space. This will help us to understand the uncertainties and 
correlations between the parameters. Furthermore, the Bayesian evidence offers a 
useful tool to select models |176l 1177] . The evidence is a normalization constant 
and is defined as: 

Z = J P(D\®,H)P(@\H)d®. (5.3) 

The evidence is independent of the parameters and therefore it is usually neglected 
in parameter estimation. However, when comparing alternative models, the evi- 
dence is the key ingredient to choose which one is better. A model which depends 
on fewer free parameters and fits better the data will have a larger evidence. A com- 
parison between two competing models Hq and Hi can be performed by comparing 
their respective posterior probabilities as follows: 

P(gi|D) P(Djgi)P(gi) Hgi) 

P(J?„|D) P(T>\H )P(H ) W P(H y [ °- ' 

where P (Hi) / P (Hq) is a priori probability ratio for models Ho and Hi and usually 
can be set to unity, Bio is called the Bayes factor and is defined as the ratio of two 
models evidences. Given the observations D, if B w > 1, model Hi is favored versus 
model Hq, and vice versa. While the \ 2 method addresses the goodness of fit, the 
Bayesian approach provides a model selection criterion. 

The main difficulty of the Bayesian approach is its very expensive computation 
cost on the calculation of the posterior distribution, and especially the Bayesian 
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evidence. To perform sufficiently fast Bayesian analysis, a publicly available pack- 
age, MultiNest |178[ 1179] . which implements a nested sampling algorithm was in- 
tegrated with GALPROP to study cosmic ray propagation models. Compared to 
traditional Markov Chain Monte Carlo (MCMC) techniques (see e.g. [T5D] ). Multi- 
Nest is highly efficient which reduces the computation time by a factor of ~ 100. 
Using a log-likelihood function, InL (0) = — l/2x 2 , MultNest directly produces the 
evidence and the posterior distribution. Once the samples of posterior distribution 
/(©) in n-dimension parameter space are generated, it is able to estimate the one- 
dimensional (ID) marginal probability P (Oj|D, H) for the parameters of interest 
9^ by integrating /(©) over all other parameters, as: 

p(e i |D,£0= /p(e|D 1 H)de 1 ...de i _ide i+ i...de n . (5.5) 

Distinguishing from the frequentist confidence interval which indicates how fre- 
quently the observed interval contains the parameters, the Bayesian credible interval 
states the degree of belief that the parameters lie inside the interval. The two-tail 
symmetric a% credible interval [ 0~ , 6^ ] can be obtained by: 

1 — off f°° 

P (©ilD, H) dQi = —- - = / P (©ilD, H) d0i. (5.6) 

Je+ 

The integration can be calculated by counting a fraction of (1— a%) /2 of the number 
of samples falling outside each side of the interval. Two-dimensional (2D) marginal 
posterior p.d.f.s are defined in a similar way. The a% credible regions are produced 
by finding out the contours in which the integration of the 2D marginal poste- 
rior density equals to a%. The best-fit parameters which maximize the likelihood 
function is also given by MultiNest as a by-product. 



5.4.1 Models and priors 

In this section, only the DR and DRC models were studied in the framework of 
Bayesian inference. This choice was made for several reasons. Firstly, without the 
need of an arbitrary break on the diffusion coefficient, the reacceleration process 
which is expected when relativistic particles scatter on magnetic turbulence, well 
describes the secondary-to-primary ratios. This can be seen from the x 2 study 
where the DR and DRC models give much smaller % 2 compared with the rather 



high values for the PD and DC models, as shown in table 5.1 Secondly, since the 
DR model has been studied widely in the literature, it is natural to choose it as 
a reference case. Thirdly, in order to understand if convection can better explain 
the data and to study the correlation between each process, the DRC model is also 
studied. 

The solar modulation potentials are included as nuisance parameters in the 
Bayesian analysis to diminish systematic effects due to uncertainties on the mod- 
ulation potentials. Since solar modulation mainly affects cosmic ray nuclei with 
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Parameters 


Prior range 


Prior type 


Propagation parameters 


D (10 28 cm 2 /s) 
S 

Va (km/s) 

AV/dz (km/s/kpc) 


[0.5, 15] 
[0.1, 1.0] 
[0, 100] 

[0, 50] 


Uniform 
Uniform 
Uniform 
Uniform 


Source parameters 


V 

N p (10- 9 cm 2 /sr/s/MeV) 


[1.7, 2.9] 
[4.54, 4.84] 


iV(2.3, 0.2) 
AT(4.69, 0.05) 


Modulation parameters 


$HEA03 (MV) 
*ACE/CRIS (MV) 
*CREAM-1 (MV) 
*Spacelab-2 (MV) 
$AMS01 (MV) 
^PAMELA (MV) 


[350, 650] 
[226, 424] 
[595, 1105] 
[280, 520] 
[315, 585] 
[350, 650] 


A/"(500, 50) 
A/"(325, 33) 
7V(850, 85) 
7V(400, 40) 
7V(450, 45) 
7V(500, 50) 



Table 5.3. Priors for propagation model parameters and solar modulation param- 
eters. The notation A/"(/i, a) is used to represent a Gaussian distribution with mean 
fi and strandard deviation a. 



energies below a few GeV, the low energy dependence of diffusion coefficient which 
was proved to be dominant over other processes may not allow any useful informa- 
tion to be extracted. Therefore, only standard models with 77 equal to unity are 
studied here. 

Based on our current knowledge, priors are specified for the free parameters 
listed in section 5.1 to restrict parameters in physically reasonable regions. The 
propagation and source parameters characterizing a model are of interest. As shown 
in table |5.3| the prior on each transport parameter is uniform to assign equal prob- 
abilities on all the possible values within the prior range. The source parameters 
are assumed to follow a Gaussian distribution with an expected mean. The solar 
modulation parameters also adopt Gaussian priors, for which the mean values are 
chosen to be the estimated value given by each experiment. As shown in equation 
|5.3| the evidence of a model depends on the priors for the parameters. If the like- 
lihoods get higher values at lower prior probability regions, the evidence will be 
suppressed. This will increase our confidence in model rejection. 



5.4.2 Results 

Identical data sets as employed in the \ 2 study are used. For all the studied 
models, the constraints on parameters and the best-fit parameters maximizing the 
likelihood are summarized in table 5.4 and the marginal posterior p.d.f.s for the 
model parameters are produced. Examples of the posterior p.d.f.s and the 68% 
and 95% credible intervals (dark and light orange, respectively) for the DR-a and 
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Figure 5.7. The ID (diagonal) and 2D (off-diagonal) marginalized post erior p.d.f.s 
of the propagation parameters for the DR-a model (as shown in table [5^4| | , The 
dark/light orange color represents the 68%/95% credible interval. The cross is the 
posterior mean, the star the best fit. 

A negative correlation between the normalization Do and spectral index 8 of the 
diffusion coefficient is seen for all the models. This can be inferred from equation 



2.10 In order to keep a roughly constant diffusion coefficient to reproduce the 
secondary-to-primary ratios, a larger Dq leads to a smaller 8. The relationship 
between the diffusion parameters (D or 8) and other propagation parameters are 
model dependent. For instance, the correlation between 8 and the Alfven velocity 
va is negative in the DR-a model but positive in the DRC-c model. If the source 
parameters are also under study, the correlation between the injection index v and 
spectral index of the diffusion coefficient 8 can be found to be negative, as well as the 
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models 



DRC-c 


DRC-b 


DRC-a 


DR-c 


DR-b 


DR-a 


Model 


2.75 ±0.13 
2.71 


3.3 ±0.4 
3.6 


3.3 ±0.5 
3.21 


3.51 ±0.05 
3.52 


3.72 ±0.16 
3.76 


6.75 ±0.15 
6.88 


_ to 

o ° 

CO 1 
s — ' to 

GO 


0.564 ±0.015 
0.563 


0.40 ± 0.04 
0.40 


0.446 ±0.021 
0.45 


0.501 ±0.008 
0.502 


0.42 ±0.04 
0.39 


0.302 ±0.007 
0.299 


O-i 


9.4 ± 1.0 
10.1 


3.3 ±2.2 
2.9 


39.7 ±2.0 
39.5 


2.7 ± 1.4 
4.1 


1.7 ± 1.2 
0.26 


35.2 ± 1.2 
36.3 


a <s 
<_ a* 


7.1 ± 1.4 
7.9 


^ H- 


33 ±8 
32.8 






1 


dV/dz 
(km/s/kpc) 


2.358 ±0.009 
2.362 


2.50 ±0.05 
2.48 


to to 


2.368 ±0.009 
2.367 


2.45 ±0.04 
2.48 


[2.3] 
[2.3] 


<5 


4.68 ±0.03 
4.70 


4.67 ±0.03 
4.65 


[4.69] 
[4.69] 


4.68 ±0.03 
4.69 


4.67 ±0.03 
4.67 


[4.69] 
[4.69] 


'o" 

3 

to 

^ ° 

g CO 

CD 

< 


-171.78 ±0.12 
2.13 


-34.94 ±0.10 
0.30 


-39.76 ±0.08 
1.70 


-184.21 ±0.12 
2.34 


-34.60 ±0.10 
0.29 


-50.48 ±0.08 
2.36 


o 

Of 

CD 
< 

a 

en 
P 
o 

CD 

>^ 
to 

O 

M-S 
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Figure 5.8. The ID (diagonal) and 2D (off-diagonal) marginalized posterior p.d.f.s 
of the propagation parameters for the DRC-c model (as shown in table \5A\ . The 
dark/light orange color represents the 68%/95% credible interval. The cross is the 
posterior mean, the star the best fit. 
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Figure 5.9. The ID marginalized posterior p.d.f.s of the solar modulation param- 
eters for the DR-a model (as shown in table [5~4| > and the relation between the solar 
modulation parameters and the model parameters. The abscissae of the 2D contours 
are propagation parameters while the ordinates are the modul ation p arameters. The 
units of the propagation parameters follow the same in figure [5~. 4. 2 1 The dark/light 
orange color represents the 68%/95% credible interval. The cross is the posterior 
mean, the star the best fit. 
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lVTndp-l 


&HEA03 

(MV) 


& ACE— cms 
(MV) 


^ CREAM— I 

(MV) 


(f) rt iin 

^ bpacelab—2 

(MV) 


& AM SOI 

(MV) 


^PAMELA 

(MV) 


DK-a 


510 ±50 
591 


342 ± 21 
367 


850 ± 90 
945 


400 ± 40 
438 


470 ± 50 
568 




1 J I \ - u 












641 ±7 
649 


DR p 
i ^ 1 1 - ( 


490 ± 50 
415 


227.2 ± 1.2 
226 


850 ± 90 
851 


400 ± 40 
373 


450 ± 50 
447 


640 ±7 
650 


DRC-a 


510 ±50 
503 


296 ± 21 
272 


850 ± 90 
884 


400 ± 40 
377 


460 ± 50 
569 




DRC-b 












622 ± 16 
644 


DRC-c 


500 ± 50 
453 


227.4 ± 1.4 
226 


850 ± 90 
862 


400 ± 40 
464 


450 ± 50 
407 


640 ±9 
649 



Table 5.5. The posterior mean with standard deviation (the first row) and the best- 
fit parameteres maximizing the likelihood (the second row) for the solar modulation 
parameters for DR and DRC models by using only B/C ratio (labelled as a), by 
using the p/p ratio plus the proton spectrum (labelled as b), and by using the B/C 
ratio, the p/p ratio and the proton spectrum (labelled as c). 



correlation between v and the normalization of the propagated proton spectrum N p . 
The former one is because a flatter injection spectrum needs to be diffuse more to 
hold the propagated cosmic ray slope that is observed. The latter one is naturally 
obtained to fit the data, i.e. a flatter injection spectrum (lower value of v) will 
deviate more from the data if it has to fit a lower value of normalization for the 
propagated proton spectrum (N p ). 

An example of the ID marginalized posterior p.d.f.s of solar modulation and 
the correlations between solar modulation parameters and the model parameters 



is shown for the DR-a model in figure 5.9 In this case, mainly §heaos and 



<&ace-cris are correlated with the model parameters. The reason is that below 
300 MeV/n the data are only from ACE-CRIS and above this energy the most ac- 
curate data sensitive to solar modulation are from HEA03. The values of §heaoz 
or $ ace-cris are positively correlated with Dq and Va, and negatively correlated 
with 6. This indicates that for a DR model, less modulation needs a smaller 6 and 
more reacceleration. 



Comparing the results in table 5.4 and table |5.1[ the constraints from the 
Bayesian analysis are consistent with the ones from the \ 2 method for the DR- 
a model and the DRC-a model. Only the DR-a model prefers the Kolmogorov 
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Data 


AlnZ 


B/C 


10.72 ±0.11 


P/p + P 


-0.35 ±0.14 


B/C+ p/p + p 


12.43 ± 0.17 



Table 5.6. The difference in log-evidence between the DRC model and the DR 
model in the Bayesian study. 



spectrum of turbulence of 5 = 1/3. For the DR-b model and DRC-b model, the 
bias on the diffusion parameters due to the dominant weight of the primary pro- 
ton flux is diminished mainly after freeing solar modulation parameters and by 
applying Gaussian priors on the source and modulation parameters. For these two 
models reacceleration is still not favored since va is estimated to be nearly zero. 
This is also the case for the DR-c model. The Bayesian results for the DR-c model 
do not deviate significantly from the ones derived from the x 2 analysis. For the 
DRC-c model, while the x 2 study converges at va — > 0, weak reacceleration with 
va = 9.4±1.0 kms -1 is allowed with Bayesian method. This is possibly because the 
best-fit $ for PAMELA is estimated as 649 MV in the Bayesian analysis, i.e. higher 
than the one used in the x 2 study (500 MV). Increased solar activity will suppress 
the overproduction of the proton flux around a few GeV caused by reacceleration. 



5.4.3 Model selection 

When selecting between two competing models Hq and H±, the evaluated Bayes 
factors indicate the strength of evidence. Empirical thresholds on the logarithm of 
the Bayes factors are InBio = AlnZ = 1.0, 2.5, 5.0, representing weak, moderate 
and strong evidence (see e.g. |177| ). 

Assuming H is the DR model and Hi is the DRC model, for each combination 
of data sets the differences in log-evidence between the two models are shown in 
table [5T6] By incorporating the proton flux with the p/p ratio, AlnZ — > (i.e. 
ln£?oi — > 0) thus no model is favored. As indicated from the x 2 study, the antiproton 
data can be sufficiently described by the PD model. The addition of reacceleration 
and convection processes is not expected to improve the description of the data. By 
using only B/C data or combining all the data together, AlnZ > 5.0. These large 
values of AlnZ strongly support the DRC model over the DR model. To explain 
all the data, the DRC model is therefore selected as the "best" one. 

The predictions of this DRC-c model for the fitted cosmic ray spectra and ra- 



tios are shown in figure 5.11 This model can not reproduce the B/C ratio below 
~ 1 GeV. As can been seen in table |5.5| and figure |5.10| the best fit value of 
^ace-cris = 226 MV in this model, as well as the posterior mean, are close to 
the lower limit value of $ ace-cris specified in its prior. This indicates that the 
value of ^ace-cris lower than 226 MV is favored and would allow to better fit 
the ACE-CRIS data. The discrepancy has already been seen in the x 2 study, which 
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*HEA03<IUV) * 4CE . CRIS (MV) *CREA.-,(MV) 




*ams»,(MV) * PM1ELA (MV) 



Figure 5.10. The ID marginalized posterior p.d.f.s of the solar modulation param- 
eters for the DRC-c model (as shown in table [5^4| , The cross is the posterior mean, 
the star the best fit. 



can only be recovered by either adding a break in the injection spectrum or con- 
sidering a low energy dependence in the diffusion coefficient. However, a break in 
the injection spectrum is difficult to explain physically. If this break is introduced, 
it either underestimates the antiproton flux (see [HI [73]), or overpredicts the B/C 
ratio at low energy (with lower value of va than that in (441 173] ) as mentioned 
in section |5.3.1[ Taking into account a nonlinear diffusion coefficient can satisfac- 
torily explain all the data, however, it has strong correlation with reacceleration 
and convection, as indicated in the \ 2 study. Including modulation potentials as 
nuisance parameters in the Bayesian analysis will increase the correlations at low 
energy. In order to concentrate the effort on understanding the processes of reaccel- 
eration and convection, the nonlinear diffusion coefficient is therefore not studied 
with the Bayesian method. Besides these two explanations, several effects may 
also be responsible for the discrepancy. The dominant weight of the proton spec- 
trum on the fitting could be the most important reason. Even reacceleration with 
va = 9.4 ±1.0 kms -1 is allowed to fit the proton spectrum as a compensation for 
higher solar modulation potential (^pamela = 640 ± 9 MV) than the value fixed 
in the x 2 study (500 MV), it is still too weak to account for the rapid decreasing of 
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boron flux below around 1 GeV. Using the force-field approximation to model the 
solar modulation may be another problem since it is too simplified and adopting 
a different $ may significantly change the estimated values of other parameters. 
The systematic inconsistencies between data sets could also cause a bias on the 
constraints since data set from different experiments were employed in this work. 



5.4.4 The electron spectrum and the positron fraction 

In addition to the cosmic ray nuclei, PAMELA also measures precisely the electron 
and positron components in the cosmic radiation. PAMELA reported the electron 
absolute spectrum between 1 GeV and 625 GeV |181| as well as the positron frac- 
tion between 1.5 GeV and 100 GeV [751 H82| . These data were not used in the 
fitting procedure described in previous sections, and therefore allow a cross check 
with the prediction of electrons and positrons from the best model. In a traditional 
scenario, cosmic ray electrons originate from SNRs and positrons are mainly sec- 
ondary production created by cosmic ray protons interacting with the ISM. The 
positron fraction measured by PAMELA increases with energy above 10 GeV, con- 
flicting with the trend predicted by secondary production. Primary sources such as 
pulsars and dark matter may give an extra contribution. Therefore, a precise and 
reliable determination of the electron flux and positron flux (or positron fraction) 
play important roles in studying primary sources. 

The primary electron injection spectrum is normalized to the PAMELA electron 
data at 70 GeV and is tuned with a power law injection index 2.72 to fit the electron 
spectrum measured by PAMELA. The electron spectrum and the positron fraction 
are calculated based on the best model and compared with the PAMELA data, as 
shown with blue color in figure |5.12| The theoretical calculation of the electron 
spectrum agrees well with the PAMELA data. It is noted that no significant break 
in the electron injection spectrum at 4 GV, as adopted in |73j which requires a 
strong reacceleration and introduces a break on the primary proton injection index 
at 100 GeV, is necessary here to describe the data. The same configuration was 
tested using the x 2 minimization method to fit the B/C ratio, the p/p ratio and 
the proton flux. As shown in figure 5.13[ using the best-fit parameters obtained 
for this model (see DR II-c in table 5.1), a broken power law with index 1.8/2.6 
below/above 4 GV is needed in the electron injection spectrum to fit the electron 
data. Otherwise an anomalous bump arises around 1 GeV, which is caused by a 
combination effect of reacceleration and energy loss. In the best model (DRC-c 



model as listed in table 5.4) under study here, the reacceleration is weak and a 
break is not required in the electron injection spectrum. 

However, the DRC-c model (with an electron injection index 2.72) predicts an 
evident lower positron fraction than the PAMELA data. The discrepancy below 
10 GeV can be due to, for example, a charge-sign dependent solar modulation. 
But above 10 GeV, the prediction which considers only secondary positrons pro- 
duced in cosmic ray spallation, shows an opposite trend of the positron fraction 
with PAMELA data. Clearly additional components are required to interpret the 
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Figure 5.11. The B/C ratio (top), the p/p ratio (midd le) a nd the proton spectrum 
(bottom) for the best model (DRC-c) as listed in table \5A\ The dark/light orange 
(or blue) color represents the 68%/95% credible interval. The blue and orange colors 
in the top figure are plotted with the posterior values of &ACE-CRIS an d &HEA03> 
respectively. The middle and bottom figures are plotted with the posterior values of 
® PAMELA- 
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Figure 5.12. The electron spectrum (top) and the positron fraction (bottom) for 
the best model (DRC-c) as listed in table [5^4] The dark/light color represents the 
68%/95% credible interval. The blue bands are calculated for the single primary 
component model with e~ injection slope=2.72. The orange bands are calculated 
for the two primary components model for which one component has an e~ injec- 
tion slope=2.68 and another one has an injection slope=2.1. The green dash 
line is plotted using the best-fit parameters derived in |73[ . Data points are the 
measurements from PAMELA [T8T1IT82] and Fermi [79]. 
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Figure 5.13. The electron spectrum for the DR II-c model as listed in table [5TT] 
The electron injection spectrum is modeled with a broken power law with index 
1.8/2.6 below/above 4 GV (red) and a single power law with index 2.6 (blue). Data 
points are the measurements of the electron flux from PAMELA [181] ■ 



raising observed by PAMELA. Since no obvious deviation from the modeled flux 
is observed in the electron data, positron production from the nearby astrophysi- 
cal sources (e.g. pulsars) and/or from the exotic sources (e.g. dark matter) may 
contribute to the extra component in positron spectrum. Both the positrons and 
the electrons are expected to be produced in equal amounts from these extra pri- 
mary sources. To reproduce both the electron spectrum and the positron fraction 
measured by PAMELA, two primary components are considered here, employing 
the same injection indices as used in |181| . i.e. a standard primary component with 
injection index 2.69 only contributing for electrons and an extra component with 
injection index 2.1 producing equal amount of electrons and positrons. The agree- 
ment between this model (plotted in orange color) and the data can be seen in figure 
|5.12| Invoking pulsars or dark matter as the extra component need more realistic 
treatment concerning the nature of the source, for example the source distribution, 
the mass of the dark matter particle, etc. 

5.4.5 What did we learn from the Bayesian study? 

From the Bayesian study, the following conclusions can be drawn: 
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• The bias caused by the dominance of very precise proton data is reduced 
in the Bayesian analysis by specifying priors on source parameters and solar 
modulation parameters, but still could possibly exist. 

• Strong reacceleration is required to explain the B/C data but still yields to 
too many protons below a few GeV, as concluded in the x 2 study. When 
including the proton spectrum to the fitting procedure, unlike va — > in 
the x 2 study, weak reacceleration is allowed in the Bayesian analysis since 
uncertainties on modulation parameters are taken into account. 

• Negative correlations are found between parameters D and 5, v and 5, as 
well as v and N p , as expected. Correlations between other parameters are 
model dependent. 

• To fit only the current PAMELA data of the p/p ratio and the proton spec- 
trum, reacceleration and convection can not be constrained. 

• When fitting only the B/C data or all the data, the DRC model is found to 
perform better than the DR model based on the Bayesian evidences. However, 
they still cannot reproduce the low energy B/C data reported by ACE-CRIS 
due to a prior limited on $ace-cris- 

• The source and propagation parameters are not able to be well constrained by 
fitting PAMELA proton and p/p data alone, i.e. errors on va and dV/dz are 
comparable to their posterior means. To simultaneously fit all the data, the 
errors on parameters range from 2% on S to 20% on dV/dz. The uncertainties 
on S are small since 6 is well constrained by the high energy data. At low 
energies, multiple physical processes shape the cosmic ray spectra and there- 
fore the estimation of the parameters characterizing low energy processes is 
hindered by their degeneracy. 

5.5 Conclusion 

Cosmic ray propagation models have been studied using both the x 2 minimization 
method and the Bayesian method. Different combinations of data sets have been 
used to test their constraining ability. Using only the antiproton and proton data 
from PAMELA does not allow to constrain propagation models and therefore B/C 
data provided by other experiments were also included in the analysis. The B/C 
data is sensitive to propagation parameters but insensitive to source parameters. 
Models with strong reacceleration or with a nonlinear diffusion coefficient can well 
reproduce the B/C data. Other models fail to reproduce the B/C ratio in the low 
energy range. In order to constrain both the propagation and the source parameters, 
a simultaneous fit on the B/C ratio, the p/p ratio and the proton spectrum have 
been performed. However, in the x 2 study, reacceleration models are disfavored in 
the simultaneous fit since reacceleration produce too many protons at low energy. 



5.5. Conclusion 
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The uncertainties due to solar modulation and the bias caused by the very pre- 
cise proton data in the fitting procedure are taken into account in the Bayesian 
analysis by specify priors on the modulation and the source parameters. The mod- 
els accounting for the low energy dependence of the diffusion coefficient were not 
studied using the Bayesian method since this effect is dominant over other low 
energy processes, i.e. reacceleration and convection. Including solar modulation 
parameters as nuisance parameters in the Bayesian analysis increases the number 
of correlated parameters and no strong constraints are expected to be obtained. 

From the estimated evidences in the Bayesian study, the diffusion model with 
weak reacceleration and convection is selected as the best one to explain all the 
data. However, a deviation between the low energy part of B/C data and the 
theoretical calculations is seen. More reacceleration or an extremely low modulation 
potential (about 100 MV) for ACE-CRIS could solve this problem. But since 
protons are much more precisely measured than other species, a small deviation 
caused by an increase of va can significantly increase the \ 2 value (or decrease 
the likelihood) and therefore will be excluded in the fitting procedure. Regarding 
using $ ~ 100 MV for the ACE-CRIS experiment, this modulation level is too low 
to characterize the solar activity although the B/C data observed by ACE-CRIS 
using in this work is taken during 1997-1998 solar minimum period. A low energy 
dependence in the diffusion coefficient may also account for the low energy B/C 
data, as studied with x 2 method. Moreover, the inconsistencies between data sets 
from different experiments could result in the disagreement between the data and 
the model prediction.lt is also worth to stress that no model studied in this work 
can account for the feature of the PAMELA proton spectrum at ~ 200 GeV, since 
a source dispersion as suggested in |168| or an acceleration model as proposed in 
|169j is not taken into account. The electron flux and the positron fraction were 
also produced by this model and compared with the PAMELA data. The model 
with a standard primary component of electrons predicts a lower positron fraction 
compared to the PAMELA data. Especially above 10 GeV the positron fraction 
measured by PAMELA increases with energy but the model decreases with energy. 
But by adding an extra primary component, both the electron and positron fraction 
can be reproduced. Nevertheless, this model can reproduce most of the observations 
except for the B/C ratio below 1 GeV. More robust constraints will be studied in 
the future by using upcoming PAMELA nuclei data. 



This study not only improves our understanding on the cosmic ray acceleration 
and propagation mechanisms, but will also provide a useful tool to study the astro- 
physical sources (e.g. pulsars) or search for exotic contributions (e.g. dark matter). 
Cosmic ray antimatter (e.g. p and e + ) are produced as secondary products in the 
ISM. They maybe also produced as the final state of dark matter annihilation and 
then propagated in the Galaxy before reaching Earth. An accurate and reliable 
estimation of the propagation will help us to discriminate whether the final state 
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includes quarks. It also allows us to constrain the mass and the cross sections of 
dark matter, since different properties on dark matter are expected to give different 
amount of contribution in the antimatter spectrum. Similarly, the diffuse gamma 
ray emissions can also be studied based on cosmic ray propagation since they are 
also expected to be produce by cosmic ray proton interactions in ISM, and will 
provide another channel of dark matter search. The cosmic ray e pairs can also 
be generated from pulsars. The fast energy loss of electrons indicates that only 
produced locally can reach Earth. Based on the determination of cosmic ray 
propagation, properties of nearby pulsars such as the injection spectrum and the 
energy conversion factor can therefore also be constrained. 



Chapter 6 

Discussion and outlook 



Studying the propagation of Galactic cosmic rays can help us understand cosmic 
ray astrophysical sources, the properties of the Galaxy including the Galactic mag- 
netic field, the interstellar medium and the nuclear interactions happening therein. 
However, since the discovery of cosmic rays one century ago, questions regarding 
the details of their acceleration and propagation mechanisms are still under de- 
bate. Different acceleration models predict different values of the injection spectral 
index. Whether and how other possible processes such as reacceleration and con- 
vection play a role in cosmic ray propagation is still not certain. The diffusion 
index S related to the spectrum of magnetic turbulence differs from 0.2 to 0.9 in 
the literature. 

Studies on cosmic ray acceleration and propagation, rely on accurate measure- 
ments on cosmic ray nuclei over a broad energy range. The determination of source 
and propagation parameters, which plays an important role for us to understand 
relevant mechanisms, are based on fitting the data of secondary-to-primary ratios 
and primary fluxes. Therefore the precision and reliability of the parameters are 
limited by the uncertainties and energy ranges of cosmic ray species measured and 
possible systemic effects existing in different experiment. To precisely measure light 
nuclei cosmic ray fluxes (e.g., protons, helium nuclei, antiprotons) and secondary-to- 
primary ratios (e.g., B/C, p/p, 2 H/ 4 He and 3 He/ 4 He) over a wide energy range, the 
satellite-borne experiment PAMELA was launched in 2006 and has been studying a 
variety of cosmic ray species for almost six years. The absence of atmospheric over- 
burden and the long live time makes the PAMELA measurements more accurate 
than those from balloon-borne experiments. These observations provide a wealth 
of opportunities to further study the acceleration and propagation mechanisms of 
Galactic cosmic rays. 

Among all the scientific objectives PAMELA are designed for, the measurement 
of cosmic ray antiprotons is one of the primary task through which not only the 
propagation mechanism but also exotic sources can be investigated. In the first 
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part of this work, the cosmic ray antiprotons were identified with the PAMELA in- 
strument over backgrounds presented by other cosmic ray components and particles 
produced by cosmic rays interacting with the experiment materials. The selection 
efficiencies of each individual detector and other correction factors concerning the 
geometrical factor, hadronic interaction losses, the live time of measurements and 
transmission through the geomagnetic field were estimated. Finally the antipro- 
ton energy spectrum and the antiproton-to-proton flux ratio were reconstructed 
between 1.5 GeV and 180 GeV. 

In the second part of this thesis, different propagation models were studied in 
this work, by using the PAMELA antiproton data including a low energy part 
down to ^100 MeV (not described in this thesis), the PAMELA proton data from 
400 MeV to 1.2 TeV, and the B/C ratio reported by previous experiments but with 
comparable energy range expected by PAMELA. The GALPROP code which solves 
the cosmic ray transport equation numerically, was employed to simulate cosmic 
ray propagation. Analyses have been performed relying on statistical methods, 
i.e. the x 2 minimization method and Bayesian inference. In previous studies, 
statistical analyses were mainly carried out by using semi-analytical models (e.g. 
most recently [1291 11301 172"] ) or were only performed for the numerical diffusion 
reacceleration model with a break in the injection spectrum index [73] . In this 
thesis different GALPROP models such as the plain diffusion (PD) model, the 
diffusion reacceleration (DR) model, the diffusion convection (DC) model and the 
diffusion reacceleration convection (DRC) model are studied for the first time based 
on statistical analyses. Models without an artificial break on the injection spectrum 
and on the diffusion coefficient are the main focus of the analyses. 

Different combinations of data sets are used to constrain models in this work. 
Using only PAMELA data is expected to minimize uncertainties due to inconsis- 
tencies between data sets, however, current PAMELA data (the p/p ratio and the 
proton flux) has been proved to be not enough to constrain propagation parameters. 
Stronger and more reliable constraints are allowed by a simultaneous fit including 
the PAMELA data as well as the B/C ratio from other experiments. 

The goodness of fit of each model was studied in the \ 2 study. Only models 
considering a low energy dependence in diffusion coefficient due to nonlinear MHD 
waves can describe simultaneously the B/C data as well as the p/p ratio and the 
proton spectrum. However, since the effect of nonlinear diffusion coefficient dom- 
inates at low energy, other processes (i.e. reacceleration and convection) are not 
possible to be studied. Models with a linear diffusion coefficient either cannot fit 
the B/C ratio below 1 GeV (PD and DC models) or generate too many protons at 
a few GeV (DR and DRC models). To reduce the uncertainties on the results due 
to solar modulation and the possible bias due to the dominance of the PAMELA 
proton spectrum in the fit the Bayesian analysis which specifies priors on the source 
parameters and the solar modulation parameters is used. The p.d.f.s of different 
parameters and the correlations between them are also able to be studied. From the 
X 2 study, the DR and DRC models can explain the B/C ratio well but cannot fit 
the data of the proton flux which is more prone to solar modulation and systematic 
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effects. These two models are studied in the Bayesian analysis by considering priors 
are specified on the solar modulation and source parameters which can reduce the 
uncertainties due to solar modulation and the possible bias caused by the dominant 
proton data in the fit. Based on the Bayesian evidences, the DRC model favors a 
Karaichnan turbulence and has been proved to be better than the DR model when 
describing all the data. The credible interval of the parameters and the fluxes 
have been shown based on the posterior samples produced in the Bayesian method. 
Only weak reacceleration and convection are allowed in the DRC model. The B/C 
ratio above 1 GeV, the proton and antiproton data can all be reproduced by the 
DRC model. The electron flux and positron fraction can also be accommodated in 
this model if an additional primary component of electrons and positrons is taken 
into account. However, the predicted B/C ratio is still higher than the data below 
1 GeV. 

Several effects, including inconsistencies between data sets, solar modulation 
and the dominance of the proton data in the fitting procedure, might influence the 
reliability of the results and result in misinterpretation: 

• The B/C ratio data used in this work are from experiments other than 
PAMELA. Systematic discrepancies may exist between data sets and result 
in a bias in the model constraints. 

• Solar modulation is a main factor affecting the spectra of cosmic ray nuclei 
and electrons below 10 GeV and tens of GeV, respectively. The simplified 
force-field approximation depending on a single parameter, the modulation 
potential is used in this work to model the solar modulation and may bias 
the results. A more realistic solar modulation should include a charge-sign 
dependency. 

• Including the proton spectrum in the fit provides constraints on the source 
parameters. However, since the proton spectrum is measured more accurate 
than data of other species, it has a dominant weight in the fit. Any systematic 
bias in the proton data may therefore significantly bias the results. 

The forthcoming secondary-to-primary ratios measured by PAMELA, including 
the B/C ratio between 100 MeV and 200 GeV, the 2 H/ 4 Hc ratio between 100 McV/n 
and 700 McV/n, and the 3 He/ 4 He ratio between 100 MeV/n and 900 MeV/n, are 
expected to allow better and more robust constraints on transport parameters. 
The B/C ratio provided by PAMELA is expected to be more precise than previous 
published data and can hopefully help clarify the longstanding issue concerning 
the value of the spectral index of the diffusion coefficient. Degeneracy between 
diffusion and other low energy processes, i.e. reacceleration and convection, is 
therefore expected to be broken. Incorporating secondary-to-primary ratios and 
primary fluxes exclusively from PAMELA, the bias due to data set inconsistencies 
and solar modulation uncertainties can be reduced. Moreover, when using data from 
a single experiment, more realistic modulation models, which require a number 
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of parameters to describe the modulation effects [43l I183| I184j , will be easier to 
treat and further decrease the uncertainties due to solar modulation. The bias 
caused by the very precise proton data is difficult to eliminate. One simple way is 
increasing the errors on the proton data but some features in the proton spectrum 
may disappear and the corresponding source information may be lost. Additionally 
the model test ability of the data will be lost. The reliability of the results can also 
be checked by seeing whether the results are consistent by fitting only the secondary- 
to-primary ratios and by fitting the secondary-to-primary ratios plus the primary 
fluxes. 

Observations in other channels, for example electrons, positrons and gamma 
rays, will give a consistency check on the acceleration and propagation models. 
PAMELA has measured the electron flux up to 625 GeV and positrons up to 
300 GeV 188 . The 7-ray diffuse emission has been recently measured by the 
Fermi-Large Area Telescope between 100 MeV and 10 GeV |189j . The statistical 
techniques applied here can also be adopted to other channels. Multi-messenger 
observations and statistical analyses on the measurements allow us to obtain com- 
plementary information on cosmic ray acceleration and propagation models. 

The AMS02 experiment, successfully installed on the International Space Sta- 
tion (ISS) in May last year, will also measure the energy spectra for a wide range 
of cosmic ray species in near furture. It has an acceptance of 0.5 m 2 sr |185j which 
is two orders of magnitude larger than that of PAMELA (21.5 cm 2 sr), and is ex- 
pected to operate on ISS for more than 10 years. The large acceptance and the 
long lifetime of AMS02 will drastically increase the statistics of the measurements 
than that have been achieved by all the previous experiments. The data with un- 
precedented accuracy may dramatically improve the model constraints and may 
potentially allow us to understand the acceleration and propagation mechanisms. 

Furthermore, based on accurate and reliable constraints on the cosmic ray prop- 
agation models, primary contributions causing from nearby pulsars or dark matter 
can be probed indirectly through anomalous antimatter components in cosmic rays 
(mainly antiprotons and positrons) and gamma rays, created during dark matter 
particle annihilation. The dark matter contribution can be extracted with respect 
to the astrophysical background of cosmic ray antimatter and gamma ray diffuse 
emission. The dark matter interpretation has been proposed to explain the positron 
excess above 10 GeV first observed by PAMELA and then confirmed by Fermi. Ad- 
ditional information on whether the dark matter annihilations result in hadronic 
final states will be given after the antiproton data at high energy available from 
AMS02. Information on the dark matter cross section and mass can also be con- 
strained by the spectra of cosmic ray antimatter and gamma rays. 

Using upcoming data to improve the constraints on acceleration and propagation 
models and to investigate properties of primary sources such as nearby pulsars or 
dark matter will be an important future task and a development of this thesis. 
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